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SUMMARY 


The characteristics method is presented for computing highly underexpanded 
free- jet boundary contours. This procedure, permitting extension of the cal- 
culations to very high pressure ratios and large distances downstream of the 
nozzle exit, has been incorporated into a computer program which facilitates 
rapid calculations of jet exhaust characteristics. Computer programs written 
in FORTRAN language are presented for the leading characteristic line, corner 
expansion fan, and the characteristic network. An experimental investigation, 
undertaken to verify the theoretical procedure, consisted of obtaining schlieren 
photographs of unheated jet exhaust boundaries from nozzles having design exit 
Mach numbers of 1.0, 2.22, and 5*0. 

Theoretical jet boundaries have been calculated and contour plots developed 
which permit construction of jet boundaries at any ratio of nozzle-exit pressure 
to ambient pressure up to 107, 10^, and 0.5 X 10^- for nozzle-exit Mach numbers 
of 1.0, 5-0 (flow angle at nozzle lip, 15 ° ), and 4.79 (flow angle at nozzle lip, 
26 . 5 °) j respectively. It was found that on logarithmic coordinates, the maximum 
boundary size varied almost linearly with pressure ratio. Reasonable agreement 
of theoretical and experimental jet plumes was obtained over the complete range 
of pressure ratios experimentally investigated. 


INTRODUCTION 


The flow of rocket exhaust gases from highly underexpanded nozzles oper- 
ating in near -vacuum conditions results in large billowing jet plumes. As a 
result, various problems arise depending on the vehicle altitude and the prox- 
imity of the rocket nozzle to adjacent surfaces, which may be either parts of 
the main vehicle or surfaces perpendicular to the jet blast, such as would be 
encountered during lunar landings and take-off. For flights of missiles within 
the atmosphere jet pluming can produce flow separation over the missile after- 
body which in turn affects the stability characteristics. (See refs. 1 to 4.) 
The radial flow direction of the exhaust gases downstream of the nozzle exit 
can exceed 90 ° relative to the nozzle thrust axis, depending upon such factors 
as nozzle-exit Mach number, ratio of nozzle-exit pressure to ambient pressure, 


and the ratio of specific heats of the exhaust products. Under such circum- 
stances, stage separation effects must be evaluated -with respect to forces and 
pressures on booster cases during the separation sequence. (See ref. 5.) 
Structural and heating problems may be introduced during orbital rendezvous 
missions in which the use of attitude control rochets for docking or trajectory 
corrections results in direct impingement of exhaust gases on parts of the 
vehicle. (See refs. 6 to 8.) 

Additional problems, created by the . interference of rocket exhaust plumes 
with communications between missiles and ground stations, have resulted in 
extensive research into the cause and effects of radio signal attenuation. A 
general discussion of several of the parameters influencing the propagation of 
telemetry signals through exhaust plumes is found in reference 9 , and a more 
limited study of a specific missile configuration is contained in reference 10. 
An extensive summary of the more recent information on signal attenuation in 
solid-propellant rocket exhaust for many of the current missiles is contained 
in reference 11, and reference 12 shows that signal attenuation also occurs 
with liquid-propellant motors. 

Jet pluming in connection with soft lunar landings and/or take-off, and the 
possibility of crater formation due to erosion of the landing area by jet 
exhaust plumes, has been the object of both experimental and theoretical 
studies. Jet impingement on hard flat surfaces has been investigated with cold 
jets (ref. 13) and with hot jets (ref. l4) . Experimental investigations of jet 
impingement on a dust-covered surface (ref. 15) indicate the possibility of 
visibility impairment and vehicle surface damage. Further theoretical studies 
of this impingement problem are contained in references 1 6 and 17. 

As a result of the many problems associated with highly underexpanded noz- 
zles, many studies have been made to obtain a clearer understanding of the 
exhaust plume behavior. A few of these studies are contained in references l8 
to 27 for both free jets and jets exhausting into a supersonic stream. Experi- 
mental confirmation of jet plume calculations by the method of characteristics 
or by empirical methods has been limited by lack of experimental data at high 
pressure ratios. 

This report presents experimental studies and theoretical studies by the 
method of characteristics of free jets exhausting into still air from sonic and 
supersonic nozzles at relatively high pressure ratios. Experimental tests, 
conducted with cold air at a total temperature of about 90° F and at stagnation 
pressures up to approximately 2,400 psia, consisted of obtaining high-speed 
schlieren photographs of the jet exhaust plume. The initial phase of this 
investigation concerning experimental tests of free jets impinging on a parallel 
flat surface at ratios of total pressure to ambient pressure up to 250,000 has 
been reported in reference 8. Theoretical confirmation of the experimental pro- 
gram required extensive modifications to the characteristics method used in 
reference 18 to obtain meaningful calculations at large distances downstream of 
the nozzle exit for high pressure ratios. In view of the continued effort 
anticipated in the study of jet plumes and the complexity in setting up the 
computation procedure, the complete program for the characteristics method is 
included in this report. 
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SYMBOLS 


dj nozzle-exit diameter 

l distance along jet axis from plane of nozzle exit to location of 

Riemann wave 

M Mach number 

Mj nozzle -exit Mach number 

Pj static pressure at nozzle exit 

• nozzle total pressure 

> J 

p^ ambient static pressure 

r perpendicular distance from jet axis, equivalent to y in appendixes 

r. radius of nozzle exit 

J 

S diameter of Riemann wave 

x distance downstream from plane of nozzle exit measured parallel to jet 

axis 

initial turning angle, measured between jet axis and tangent to jet 
boundary at the nozzle lip 

7 ratio of specific heats 

0 angle of flow, measured between direction of flow and jet center line 

0$ nozzle half -angle , or angle of flow at nozzle lip, measured between 

direction of flow and jet center line 

p Mach angle 

v Prandtl-Meyer angle (angle through which a supersonic stream is turned 

to expand to a Mach number greater than l) 

Vjy Prandtl-Meyer expansion angle from sonic velocity to nozzle -exit Mach 

number 

Prandtl-Meyer expansion angle from sonic velocity to jet-boundary Mach 
number 

Subscript: 

max maximum 
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APPARATUS FOR EXPERIMENTAL STUDIES 


Test Setup and Procedure 

The principal part of the experimental investigation to obtain schlieren 
observations of free -jet flow fields was conducted in the 4l-foot-diameter 
vacuum sphere at the Langley Research Center with the test setup as shown in 
figure l(a). Air from a tank farm pressurized to approximately 2,400 psia and 
90° F, and of sufficient volume to maintain essentially a constant nozzle stag- 
nation pressure during a test run, was supplied to the nozzles through a 
3/4-inch-diameter supply pipe. An electrically operated solenoid valve located 
just upstream of the nozzle permitted a rapid start sequence. Over a total 
testing time of approximately 30 seconds the ratio of total pressure to ambient 
pressure was reduced from a value of approximately 250,000 down to about 20,000. 

Nozzle stagnation pressures were measured by a 3^000-psia pressure trans- 
ducer located between the electrically operated solenoid value and the nozzle 
inlet bell. (See insert in fig. 1(a).) Initial pretest sphere pressures of 
approximately 0.4 mm Hg (0.0077 psia) were attained with vacuum pumps and meas- 
ured with a McLeod gage. The change in sphere ambient pressure with time was 
measured by a pressure transducer calibrated from 0.005 to 0.100 psia. All 
pressure measurements were continually recorded on oscillograph paper. With the 
test nozzle in operation, the sphere pressure increased linearly with time; 
therefore, the ratio of total pressure to ambient pressure decreased inversely 
with time. 

High-speed l6-mm motion pictures of the exhaust plume were obtained from a 
double-pass schlieren system. A 25-inch-diameter parabolic mirror was mounted 
about 4 feet behind the nozzle (as shown in fig. 1(a)), and both the camera and 
schlieren light source were located on the equator of the vacuum sphere in the 
monitoring room. The camera was synchronized with the flashing light source and 
operated at about 750 frames per second. Timing marks recorded on the edge of 
the film from a 60-cycle flashing light source permitted a correlation of the 
motion pictures with a plot of pressure ratio against time. 

Supplementary experimental data were obtained in the low-pressure-ratio 
test facility shown in figure 1(b). This test setup permitted operation at a 
constant pressure ratio and adjustment of the single-frame schlieren system for 
maximum sensitivity. The larger nozzles tested in this facility improved the 
quality of the schlieren photographs of the free- jet flow fields. Supply air 
from a large -capacity tank farm was available at total pressures up to 400 psia 
and total temperatures of about 90° F. Stagnation pressure was measured by a 
total-pressure tube located at the outlet of a small settling chamber as shown 
in the figure. With a nozzle in operation, values of ambient pressure as low 
as 0.5 psia could be held constant over long periods of time by using a series 
of compressors as vacuum pumps. Ambient pressure in this facility was measured 
with a mercury manometer. 
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Test Nozzles 


Sketches of the two nozzles used for the high-pressure-ratio tests in the 
4l-foot vacuum facility are shown in figure 2(a). The converging nozzle 
(Mj = 1.0) had an exit diameter of 0.125 inch, and the converging-diverging 

conical nozzle, having a nominal design Mach number of 5-0, based on inviscid 
flow, had an exit diameter of 0.625 inch, an area ratio of 25, and a half -angle 
of 15°. For some tests conducted with the nominal M = 5*0 nozzle, a static- 
pressure orifice installed in the diverging wall just upstream of the exit indi- 
cated an actual exit Mach number of 4.79* 

The initial phase of this investigation (ref. 8) reported that measure- 
ments of the initial turning angle of the flow at the nozzle exit a^, obtained 

from a series of photographic enlargements of the Mach 5 nozzle, indicated a 
substantially larger value of expansion half -angle 0^ than the inviscid design 

value. This comparison was performed by using an expression, on = V 1 ~ V N + 
(see, for example, ref. 19). With and determined from the estimated 

exit Mach number of 4.79 and the experimental value' of the ratio of ambient 
pressure to total pressure, respectively, the effective nozzle half -angle 0^ 
was computed to be 26. 5° * or about 11.5° greater than the inviscid design value. 
The difference between the measured and calculated value of a^, when correlated 

with a Reynolds number of about 1J x 10^ based on nozzle-exit conditions and 
diffuser conical length, compares favorably with the results of reference 5- 
The larger effective nozzle half -angle is believed to be due to a combination 
of factors, including small-scale model size, sharp corners at the nozzle 
throat, boundary-layer growth, and a thinning of the boundary layer at the noz- 
zle exit. It is evident that more research is needed in this area; and until 
experimental results indicate otherwise, it is believed that 0^ will corre- 
spond to the geometric angle for large-scale models. 

The two large-scale models tested in the low-pressure -ratio facility are 
shown in figure 2(b). Both nozzles had smooth entrances and minimum throat 
diameters of 0.707 inch. The Mj = 2.22 contoured nozzle was designed by the 

method of characteristics for parallel flow at the exit. Nozzle coordinates 
are tabulated on the figure. 


RESULTS AND DISCUSSION 


Basic Considerations 

Jet exh aust description .- Prior to a discussion of theoretical and experi- 
mental results, a few of the more basic fundamentals of jet exhaust plumes will 
be considered. Rocket nozzles are usually classified either as conical or con- 
toured. The exhaust gases leave a conical nozzle in a radial direction whereas 
the contoured nozzle is designed for exhaust flow parallel with the axial center 
line. Regardless of the type of nozzle used the structure of the exhaust jets 
produced by operating these nozzles in a highly underexpanded condition is ^ 
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basically similar , as shown by the illustration in figure 3 • Some of the more 
important features of a typical exhaust plume consist of the leading character- 
istic line, internal shock, jet boundary, and normal shock. 

The leading characteristic, or Mach line, as determined from the nozzle- 
exit Mach number, is a straight line for the contoured nozzle and forms the 
familiar Mach cone. All flow conditions along this Mach line and within the 
Mach cone, including M, 0, p, and v, are constant. For the conical nozzle 
the leading characteristic line is curved and flow conditions are different for 
each point downstream of the nozzle exit. A more detailed discussion of the 
leading characteristic line is contained in appendix A. For both nozzles, how- 
ever, when operating at or above the design pressure ratio the leading charac- 
teristic serves as a dividing line between two different flow regions. The 
internal region (Mach cone or radial flow zone) bounded by this line consists 
of a zone of zero external influence in which flow properties are determined 
strictly from nozzle design criteria. In the region external of this line and 
bordered by the jet boundary, the ratio of stagnation pressure to ambient pres- 
sure is the controlling factor with respect to size and shape of the jet plume. 

In the region bounded by the leading characteristic line, the axial center 
line, the normal shock, and the internal shock, the controlling factors are the 
corner expansion fan originating at the nozzle lip and the overall three- 
dimensional nature of the flow. The characteristic net comprising this region 
is made up entirely of an expansion network in which the flow is isentropic. 
Mach numbers within this area are therefore continually increasing, both with 
distance downstream of the nozzle exit and also in a direction perpendicular to 
the axial center line out to the position of the internal shock. 

The internal shock is formed by the coalescence of infinitesimally weak 
compression waves reflected from the jet boundary. At pressure ratios slightly 
above the nozzle design value of 1.0 the internal-shock system forms the famil- 
iar diamond pattern. As pressure ratio increases a normal shock forms and the 
internal shock approaches the jet boundary. The internal shock originates near 
the nozzle exit and increases in strength with distance downstream of the noz- 
zle. In traversing the internal shock in a direction perpendicular to the noz- 
zle axial center line, a sharp Mach number discontinuity exists followed by a 
gradual decrease to the final boundary value. The jet boundary, as the name 
implies, defines the external shape of the exhaust plume and is the separation 
line between ambient conditions and flow from the exhaust nozzle. Mixing at the 
jet boundary is ignored in this presentation. 

Factors affecting exhaust plume .- The primary parameters which influence 
exhaust plume size and shape are Mj, 0jj, y, Pj, and p cx) . Of these param- 
eters Mj , 0jj, and pj are constant and dependent only upon nozzle design 

parameters and chamber stagnation pressure. In real-gas flow the ratio of 
specific heats y varies to some degree with the reduction in exhaust gas tem- 
perature produced in the expansion process, and the state of equilibrium of the 
flow. In order to prevent further complications of the calculation procedure, 
frozen equilibrium has been assumed; therefore, y remains constant. The 
remaining parameter is ambient pressure p^ which varies with altitude and is 

expressed as nondimen sional pressure ratio, Pj/p M . 
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Highly underexpanded exhaust nozzles have an excess of pressure pj at 

the nozzle exit, and, as a result, the exhaust gases leaving the nozzle exit 
must expand until the pressure is equal to the surrounding or ambient pressure. 
At near vacuum values of ambient pressure, even moderate nozzle -chamber pres- 
sures produce very large values of p^p^* It is at these high pressure ratios 

that concern arises as to the maximum angle through which the exhaust gases 
expand and the overall size and shape of the plume. 

Initial turning angle .- The method of calculating the initial turning 
angle at the nozzle lip is shown in figure 4 along with a sketch of the initial 
section of the jet boundary. Also listed are the controlling parameters pre- 
viously discussed. Since the circumference of the nozzle exit can be considered 
to be made up of a large number of straight line segments, the flow just at the 
nozzle lip can be classified as two dimensional. Furthermore the absence of 
shocks at the nozzle exit indicates an isentropic expansion of the flow to the 
ambient pressure. The flow angle through which the flow will expand 

is the difference between the Prandtl-Meyer expansion angles 

corresponding to the jet-boundary Mach number and the nozzle-exit Mach number. 
Boundary Mach number, the only unknown in determining v^, depends upon the 

altitude, with its corresponding ambient pressure p ro , and the nozzle total 
pressure Pt,j* Inclusion of the flow direction as imparted by the nozzle half- 
angle completes the equation for finding the initial turning angle: 


a N “ V 1 " V N + 0 N 


An indication of the magnitude of the initial turning angle and the varia- 
tion with pressure ratio is shown in figure 5 for two different nozzles. For a 
perfect gas (7 = 1.4) the maximum Prandtl-Meyer angle through which the flow 
might expand at an infinite pressure ratio is shown in reference 28 to be 
130.45°. The expansion process may be either internal in which part of the 
expansion occurs within the nozzle proper to attain a given exit Mach number, 
or external in which the process is concluded by expanding to ambient conditions 
after the gas has left the nozzle. In order to achieve the maximum turning 
angle (130.45°) all of the expansion is required to occur external of the noz- 
zle as is shown for the Mj = 1.0 case in figure 5 } however, this situation is 

not desirable, since external expansion is synonymous with loss in thrust per- 
formance. From the standpoint of increasing propulsive efficiency, it is 
desirable to increase the Mach number of the nozzle in order to achieve as much 
internal expansion as possible. Structural and weight limitations, however, 
place a limit on the nozzle design Mach number, which in turn precludes having 
all internal expansion, particularly for high-altitude operations. As a result 
of these limitations, nozzles are usually designed for a moderate supersonic 
Mach number such as the Mj = 5-0 case shown in figure 5 . The maximum turning 

angle for this particular design Mach number and nozzle half -angle is only 
about 69 ° • 
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If the Mj = 5.0 nozzle of figure 5 were contoured for axial flow at the 
exit, 0N = 0°, the variation of initial turning angle with pressure ratio would 
he a curve originating at = 0°, similar to the existing curve, and located 
15° lower at any selected pressure ratio. 

It would appear from the previous discussion that restrictions of the 
turning angles to moderate values might not he too difficult. When real-gas 
flows are encountered the initial turning angle may he much higher than those 
indicated in figure 5* As an example, initial turning angles with all external 
expansion reach maximum values as high as 243° for 7 = 1.15, and 203° for 
7 = 1.20 (see ref. 29), in comparison with 130.45° for 7 = 1.4. Although some 
of the expansion process will he internal, as has previously been shown, these 
numbers serve to illustrate that initial turning angles approaching and even 
exceeding 90°, depending upon the value of 7, should not he unexpected. 


Theoretical Results 

Method of calculations .- The method of characteristics used in calculating 
free-jet boundaries, by using three-dimensional irrotational equations of flow 
and the "foldback" procedure, has been well established in reference l8. How- 
ever, certain limitations arise in using the foldback method as the pressure 
ratio increases. These limitations restrict the range- of computations to a 
region close to the nozzle exit. In order to obtain complete contours of the 
jet boundary out to the maximum diameter of the jet plume at the higher pres- 
sure ratios, it became necessary to eliminate the foldback portion of the cal- 
culations and develop a procedure in which the outer region of the character- 
istic net would be developed similar to that of the actual flow. A description 
of the calculation procedure including equations and a complete FORTRAN (FORmula 
TRANslation) writeup is included in the appendixes. Consequently, only a brief 
description of the revised characteristic network obtained by eliminating the 
foldback is included here. 

A segment of a typical characteristic network is shown in figure 6 and 
illustrates points of primary interest in the revised calculations procedure. 

The nozzle-exit center line is located at the top of the figure and flow from 
the nozzle is from left to right. The initial turning angle at the nozzle exit 
CL-§ is in excess of 90° and the greatly expanded abscissa scale permits improved 
visualization of the characteristic network. Three primary types of lines com- 
prise the characteristic network and are identified in the figure as A, B, and 
C lines. Type A lines originate at the nozzle lip and comprise the corner 
expansion fan, whereas type B lines originate either on the leading character- 
istic or on the diametrically opposite lip. The type C lines are referred to 
as reflected rays due primarily to the fact that the B lines strike the boundary 
and then appear to reflect as though from a physical surface. The coalescence 
of these reflected rays, which are actually very weak compression waves, forms 
the internal shock previously discussed. As more of these waves coalesce, the 
strength of the internal shock increases. 

Additional intersections of the comer expansion rays, A lines, with the 
internal shock, in addition to the reflected ray intersections, result in 
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obvious complications to the calculations in the region of the internal shock. 

In the foldback procedure the corner expansion rays are allowed to continue on 
past the location of the internal shock and hence no solution is required for 
intersecting lines; however , this situation does not simulate actual flow con- 
ditions. A reasonable approximation of actual flow conditions is to combine 
these rays at the point of intersection. In three-dimensional calculations rays 
are generally referred to as right-or-left-running characteristic lines (or mem- 
bers of opposite families), and by using the known flow conditions on two of 
these lines it is possible to compute a third point to complete the net. From 
the sketch in figure 6 it is observed that the intersections which must be 
handled are composed of members of the same family type of lines and must there- 
fore be handled as such. Development of appropriate equations for the solution 
of flow conditions at these intersections is presented in appendix C. 

Initially, in setting up the program for machine calculations, only two 
intersections were considered possible in between any two consecutive B lines. 

It is clear, however, that the spacing of these two consecutive characteristic 
lines determines to some extent the number of intersections that might occur. 
When more than two intersections occurred the reflected rays, which would pro- 
duce the additional intersections, were skipped and the calculations proceeded 
as though these rays were nonexistent. Omitting these points amounted to 
increasing the size of the characteristic net at this particular location with 
a subsequent small decrease in accuracy. Some of these dropped points are shown 
in figure 6 as reflected rays that are terminated before intersecting the inter- 
nal shock. The FORTRAN program as presented in appendix C has been revised, 
however, to permit any number of intersections to occur between two consecutive 
B lines. 

Jet boundary coordinates for Mj = 1 as obtained from the revised cal- 
culations are compared in figure 7 with results obtained by using the foldback 
method of reference l8. Coordinates are presented in a non dimensional form, 
xj rj and r/rj, for pressure ratios Pj^/P^ of 20 and 45,000. Agreement 

between the two methods is good at VjJVoo = 20, but a decided deviation between 
the two curves occurs at large values of T j T ^ for Pj^Poo = 45,000. 

. The oretical jet boundaries .- Results of the calculations of jet plume 
boundaries are shown in figures 8(a) to 8(g) for the following nozzles: 

Mj = 1*0, % = 0°; Mj = 5.0, e N = 15°; Mj = 4.79, 0N = 26.5°; and Mj = 2.22, 

= 0°. Two sets of curves are shown for some of the test nozzles - one set, 
such as in figure 8(a), at an expanded scale to show greater detail of the 
boundaries near the nozzle axis, and the second set showing the complete bound- 
aries (fig. 8(b)). Symbols on some of the curves show point spacings on the 
boundaries but are omitted for the lower pressure ratios or where crowding of 
the curves would occur. Figure 8(b) indicates that the maximum radius of the 
jet exhaust plume will be about 450 nozzle radii for a pressure ratio of 45,000, 
with the location of the maximum radius occurring approximately 1,500 nozzle 
radii downstream. As a result of limitations of the initial program setup, 
boundaries for the two highest pressure ratios were terminated at a radius of 
approximately 550 nozzle radii. As additional boundaries were computed and 
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insight gained in the intricacies involved in the calculation procedure, these 
limitations have gradually been removed. 

Similar boundaries computed for a Mj = 5-0 nozzle with 0jj = 15° are 

shown in figures 8(c) and 8(d). In the calculations for this nozzle as well as 
all others, pressure ratios Pj/Poo were chosen from preselected values of 

initial turning angle ajj on the basis of maintaining a relatively uniform 
increase in a,jj. The maximum angle through which the gas could expand for the 
Mj = 5*0 nozzle was only about one-half of that for the previous Mj = 1.0 

nozzle; therefore, fewer boundaries were computed. At a pressure ratio of 
8,143 the boundary reached a maximum radius of 225 nozzle radii at a distance 
of 1,050 radii downstream of the nozzle exit. 

Boundary coordinates for a Mj = 4.79 nozzle with a half -angle 0^ of 
26 . 5 ° are shown in figures 8(e) and 8(f) for pressure ratios ranging from about 
29 up to 2,926. The maximum computed boundary at Pj/p^ = 2,926 corresponded 

to r jr j = 188 

The two boundaries in figure 8(g) were computed for comparison with experi- 
mental results from a Mj = 2.22 nozzle having parallel exit flow at the two 

specific pressure ratios shown. Also given are the internal shock coordinates 
as obtained from the computer program. These shocks appear to originate at some 
distance downstream of the nozzle exit; however, if an infinitely small point 
spacing is selected, the shock would approach the nozzle exit. 

Maximum jet boundary coordinates .- Rapid determination of the size of the 
exhaust jet at the maxi mum diameter of the exhaust plume may be obtained from 
figures 9( a ) "to 9( c ) for three of the nozzles for which boundaries are pre- 
sented. The maximum boundary coordinates, (x/r^) and (v/rA , are shown 

\ ' ‘vmax \ ' °/max 

as a function of pressure ratio on logarithmic coordinates. From this figure 
the maximum size of the jet boundary may be obtained for any intermediate pres- 
sure ratio over the range of computed values. The nearly linear variation on 
log-log paper of the maximum boundary size with pressure ratio permits extrap- 
olation to much higher pressure ratios with an expectancy of good accuracy. • 


at x^r j = 720 . 


Boundary contours from cross plots .- Construction of additional jet 
boundary contours at any intermediate pressure ratio, other than those specifi- 
cally computed, may be obtained from the plots of figures 10(a) to 10(c) for 
nozzles with exit Mach numbers and half -angles as follows: Mj = 1.0, 0jj = 0°; 

Mj = 5-0, 0jj = 15°; and Mj = 4.79, 0 N = 26. 5° • These plots were obtained by 


cross plotting figures 8(a) to 8(g) at constant values of 


x / r J- 


The pressure 


ratios corresponding to the computed boundaries are indicated by the ticks at 
the bottom of each figure. It is noted that all curves of constant x/rj fair 

asymptotically into a single curve which has been shown in figure 9 as the var- 


iation of 


(7 r j). 


with pressure ratio. 


10 



The only calculated boundaries in which the jet expanded in excess of 90° 
were obtained for the Mj = 1.0 nozzles, and cross plots for this condition 

are shown in figure 10(a). The initial turning angle first reached 90° at a 
pressure ratio of 1 , 852 , and this is indicated by the intersection of the line 
for = 0 with the pressure-ratio axis. At higher pressure ratios double 

values of the boundary coordinate 


r / r j 


exist for all negative values of 


x^r j. The locus of maximum absolute negative values of x / r j (in terms of 
r^/rj and pressure ratio), as defined by a line passing through the minimum 
pressure-ratio point for each different negative value of x / r j; is shown to be 
an essentially linear variation on log-log paper. At pressure ratios p./p 


on the order of 10^, values of x / r j in excess of -15 nozzle radii exist 


at values of 


r / r j 




of about 150 nozzle radii. 


Experimental Results 

Comparisons of experimen tal and theoretical results .- One of the primary 
objectives of this investigation, in addition to computing the jet boundaries, 
was to provide experimental confirmation of the calculation procedure. Experi- 
mental results therefore consist of selected schlieren photographs of the free- 
jet flow field obtained at constant flow conditions and also individual enlarge- 
ments from high-speed schlieren motion pictures in which flow conditions were 
continually varying. Photographs of the jet plume are shown in figures 11, 12, 
and 15 for the Mj = 1.0, Mj = 2.22, and nominal Mj = 5*0 nozzles, respec- 
tively. The computed jet boundaries, indicated by dots on the photographs, were 
obtained from the cross plots of figure 10 at the jet pressure ratios corre- 
sponding to those of the photographs. The series of photographs for the 
Mj = 1.0 nozzle in figures 11(a) to ll(h) and the Mj = 2.22 nozzle shown in 

figures 12(a) and 12(b) was taken in a test facility in which operating condi- 
tions could be held constant; consequently, the schlieren sensitivity could be 
adjusted for best results. Conversely, results obtained from the high-speed 
motion pictures required a presetting of the schlieren sensitivity, which did 
not always prove to be the correct setting. As a point of reference the grid 
wires shown in figures 11(a) to ll(h) and figures 12(a) and 12(b) were spaced 
for 5-inch squares; for the small-scale models, figures ll(i) to ll(p) and fig- 
ure 13 , a 6-inch spacing of the grid wires was used. 

Comparisons of computed boundaries with experimental results for the 
Mj = 1.0 nozzle in the pres sure -ratio range of pj^/p^ = 11-5 to Pj/p^ = 19^-2 

(figs. 11(a) to 11(h)) show that excellent agreement of the boundaries is 
achieved. In each of these photographs the last point on the theoretical curve 
corresponds to the location of the maximum diameter of the jet plume as obtained 
in the calculations. Downstream of the maximum plume diameter the main flow 
field appears to remain essentially constant in size; however, there is a con- 
siderable increase in the amount of turbulence along the edge of the jet. 
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For all pressure ratios for which photographs of the flow field are shown, 
the pressure ratios are sufficiently high that the familiar shock diamonds, or 
periodic chain-like jet structures, have heen reduced to a single lobe with the 
accompanying normal shock. Downstream of the normal shock the jet plume appears 
to consist primarily of turbulent subsonic flow. That the flow is subsonic is 
evident by the lack of any additional shock structure or lobes typical of those 
generally associated with free jets at pressure ratios near the nozzle design 
value. The schlieren sensitivity and viewing distance downstream of the initial 
normal shock are believed to be sufficient for detection of other lobes or shock 
structure if they were present in the flow. 

Photographic enlargements of individual frames from the motion -picture 
film, showing the free-jet discharge from the small-scale model, are shown in 
figures ll(i) to ll(p) for the Mj = 1.0 nozzle at pressure ratios ranging 

from VjjVoo = 9A90 to Pj/Poo = 102, l 8 o. Again, good agreement exists between 

the computed and experimental jet boundaries. The distance over which the 
clearly defined boundary exists appears to be a function of pressure ratio; the 
higher the pressure ratio the greater the length of boundary definition. As 
shown in figure ll(p) the boundary definition extends out to the approximate 
location of the second grid wire, 12 inches or about 96 nozzle-exit diameters, 

where it begins to form a general fan shape and becomes less clearly defined 

with increasing distance from the nozzle. At the outer edge of the photograph 

the boundary is completely obscured. The fan region appears very unstable in 

the motion pictures with rather violent in-and-out fluctuations which increase 
with distance away from the nozzle. These violent fluctuations in the neighbor- 
hood of the jet boundary were also observed and reported in reference 10 , 
described as a turbulent mixing region, in which motion pictures were made of 
the exhaust region from the rocket motor of a second-stage flight vehicle in 
the presence of an external flow field. Motion pictures of the current tests 
show considerable turbulence with the computed boundary defining the outermost 
edge of the fan region. 

An additional point of interest concerns the location of the normal shock 
with respect to the location of the maximum jet plume diameter. At low pressure 
ratios the normal shock is located downstream of the maximum jet plume diameter, 
whereas at the higher pressure ratios the shock location occurs ahead of the 
maximum diameter. A comparison of figures 11(a) and ll(i) indicates these 
phenomena . 

Schlieren photographs of free-jet flow from a Mj = 2.22 contoured nozzle 

are shown in figures 12 (a) and 12 (b) in which both the computed boundary and 
internal shock are indicated. A limiting minimum value of ambient pressure 
restricted the maximum pressure ratio to a relatively low value of Pjy/Poo = 37-7* 

Good agreement is indicated for both the jet boundary and the location of the 
internal shock. 

Theoretical and experimental jet plume boundaries for the nominal Mj = 5.0 

nozzle are shown in figures 13(a) to 13 (e). It should be emphasized here that 
the theoretical calculations were based on the actual exit Mach number of 4.79 
and nozzle half -angle of 26 . 5 °, which are both considerably different from the 


12 



inviscid design value. The good agreement between theoretical and experimental 
boundaries indicates that the method of obtaining the nozzle half -angle in the 
presence of a boundary layer of unknown thickness is an acceptable procedure 
for small-scale models. 

Further details of the flow phenomena associated with exhaust jets, which 
were not readily apparent in the individual motion-picture frames (figs. 13 (a) 
to 13(e)), were determined from the schlieren motion pictures. Boundary defini- 
tion at Pj/Poo = figure 13 (a), extended about halfway between the end of 

the nozzle and the upper edge of the photograph ^ x / r j = 29 j - Downstream of 

this approximate location severe turbulence exists on, and obliterates, the jet 
boundary. As pressure ratio increases the boundary definition extends further 
downstream of the nozzle exit and appears to delay the onset of the turbulence. 
The edge of the jet is again subject to violent in-and-out fluctuations, some- 
what less severe than those described for the Mj = 1.0 jet. Vortices appear 

to form on the jet boundary at the approximate location of the fadeout in bound- 
ary definition with substantial growth occurring with increasing distance down- 
stream of the nozzle exit. The theoretical boundary defines very closely the 
outer edge of both the visible boundary near the nozzle exit and also the tur- 
bulent region that exists further downstream. 

Comparisons of characteristic and empirical solutions .- In view of the 
obvious complications and time involved in the characteristic calculations, an 
empirical method of computing jet boundaries proves to be extremely beneficial. 
Such a method was developed in reference 19 and permits both rapid and satis- 
factory calculations of the initial portion of jet boundaries for a wide range 
of nozzle Mach numbers, specific-heat ratios, and pressure ratios. A compari- 
son of boundaries computed by the method of characteristics and by the empirical 
method of reference 19 is shown in figures l4(a) and l4(b) for the = 1.0 

and Mj = 5*0 nozzles, respectively, to determine to what extent the empirical 
solutions are in agreement. 

Results are presented in figure l4(a) for the Mj = 1.0 nozzle at pres- 
sure ratios ranging from Pj^/Poo = 45 # to pj/p^ = 45,000. In general, the 

empirical boundaries fall slightly below the characteristic curves near the 
nozzle exit but finally cross over and continually deviate in such a manner as 
to produce boundaries larger than those indicated for the characteristic solu- 
tion. Satisfactory comparisons of the boundaries are obtained for values of 
r/rj up to about 5 at the lowest pressure ratio and up to about 10 for 

Pj^Poo = 45,000. Thus, when interest in the jet boundary is restricted to a 

region relatively close to the nozzle exit, the empirical solution provides a 
rapid and relatively accurate method for approximating the boundary. Similar 
results are shown in figure 14(b) for the Mj = 5-0 nozzle at pressure ratios 

from Pj^Poo = 20 to Pj^Poo = 2,000. Correlation of the two methods of calcu- 
lations is again good with the range of agreement extending over a somewhat 
greater range than the Mach 1.0 results. 
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Location and diameter of R iemann wave (normal shock).- The effects of jet 
pressure ratio on the nondimensional distance along the jet axis from the plane 
.of the nozzle exit to the location of the Riemann wave ^ l /djj are shown in fig- 
ure 15(a) for the Mj = 1.0 nozzle. Also shown for comparison are the results 

obtained from figure 8(a) of reference l8. With both the nondimensional dis- 
tance l /dj and pressure ratio Pj^Pco plotted on log coordinates the result 

is a near linear variation over the range of pressure ratios for which the 
shock position could be determined. 

The schlieren photographs of figures 11(a) to ll(p) show that the thick- 
ness of the shock at times may be as much as several nozzle-exit diameters. 
However, this is probably a three-dimensional effect which produces an illusion 
of considerable shock thickness. The values shown in figure 15 are the average 
distances. It is of interest to note also that the motion pictures show a high- 
frequency pulsation of the Riemann wave along the nozzle axial center line with 
an amplitude roughly equivalent to the apparent shock thickness shown in the 
photographs. Instability at the outer edges of the Riemann wave diameter is 
even more pronounced with fluctuations in location occurring first on one side 
of the center line and then on the other. These fluctuations appear to be 
dependent on local pressure and become increasingly severe at the extremely low 
pressures and high Mach numbers within the exhaust plume. 

Figure 15(b) presents results of the normal-shock location for a contoured 
nozzle with an exit Mach number of 2.22. Figure 15(c) presents previously 
unpublished data from a contoured Mj = 2.65 nozzle. Again on log-log paper, 

linear variations of the shock location with pressure ratio are shown, which 
indicate that extrapolations may be made to higher pressure ratios with a prob- 
ability of fair accuracy. However, in these figures it is shown that ambient 
pressure has a definite effect on the shock location. At a given pressure ratio 
the lower the ambient pressure the less the distance from the nozzle exit to the 
location of the normal shock. Although the effect is small at the lower pres- 
sure ratios, in terms of nozzle diameters, slightly larger effects are indicated 
for the higher pressure ratios. 

The effect of jet pressure ratio upon the nondimensional diameter of the 
normal shock (s/djj is shown in figure 16 for the sonic nozzle. The current 

data span the limited data included from figure 9( a ) of reference 18 and extend 
the pressure-ratio range for which the shock diameter could be determined to 
Pj/Po, of about 50,000. In general, agreement of both sets of data is very 

good. Some uncertainty exists in measuring the shock diameter at the highest 
pressure ratios (see figs. ll(i) to ll(p)), resulting in the scatter of data 
points at the high-pressure-ratio end of the curve. The maximum scatter is on 
the order of approximately 15 nozzle diameters. The extension of data provided 
in the current tests indicates, for pressure ratios in excess of Pj/Poo 10; 

that the variation in shock diameter is essentially a linear variation with 
pressure ratio on log-log paper. At a pressure ratio of 50,000 the shock has a 
diameter of approximately 100 nozzle diameters, which is equivalent to roughly 
50 percent of the maximum jet plume diameter (see fig. 9( a ) for jet plume size). 
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SUMMARY OF RESULTS 


An experimental and theoretical investigation of highly under expanded free- 
jet boundaries has been conducted. The experimental program, initiated to 
verify the theoretical calculations of the boundaries and internal-shock posi- 
tions, was conducted with unheated air and consisted of obtaining schlieren 
photographs of exhaust jet boundaries from nozzles having design exit Mach num- 
bers of 1 . 0 , 2 . 22 , and 5 . 

Extensive modifications to an existing characteristics -method computer 
program have been made which permit extension of the calculations of free-jet 
boundaries to high pressure ratios and large distances downstream of the nozzle 
exit. The calculations procedure has been incorporated into a computer program 
which facilitates rapid calculations of jet exhaust characteristics. A FORTRAN 
program for use on an IBM 709^ electronic data processing system is presented. 

Based on the calculations using the method of characteristics and the 
experimental tests, the following results were obtained: 

1. Theoretical jet boundaries have been calculated and contour plots devel- 
oped which permit construction of jet boundaries for the following nozzles: 

Mach 1.0, Mach 5.0 (wall angle 15°), and Mach 4.79 (wall angle 26.5°) • Partial 

76 

or complete jet boundaries may be constructed at pressure ratios up to 10 , 10 , 
and 0.5 X hA for the three nozzles, respectively. 

2. For a given nozzle configuration it was found that the maximum jet 
boundary size varied almost linearly on logarithmic coordinates with the ratio 
of nozzle-exit pressure to ambient pressure. For the Mach 1.0, 5-0, and 4.79 
nozzles the calculated maximum radius of the jet boundary was 450, 225 , and l 88 
nozzle radii at a distance of 1 , 300 , 1 , 050 , and 720 nozzle radii downstream of 
the nozzle exit, respectively, at corresponding pressure ratios of 45,000, 

8,143, and 2,926. Experimental results indicated an almost constant jet plume 
size downstream of the maximum plume diameter. 

3 . Comparisons of theoretical and experimental jet boundaries show good 
agreement over the complete range of pressure ratios investigated. Limited 
comparisons of the theoretical and experimental internal- shock location also 
showed very good correlation. 

4. Experimental results at low pressure ratios show the Riemann wave (nor- 
mal shock) to be positioned downstream of the maximum plume diameter, whereas 
at the higher pressure ratios the shock location occurs upstream of the maximum 
diameter. At the higher pressure ratios the familiar shock diamond structure 
reduces to a single lobe, or normal shock. The location of the shock down- 
stream of the nozzle exit and also the maximum diameter of the shock have a 
near linear relationship on logarithmic coordinates with the ratio of nozzle- 
exit pressure to ambient pressure. 

Langley Research Center, 

National Aeronautics and Space Administration, 

Langley Station, Hampton, Va., February 26, 1964. 
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APPENDIX A 


DETERMINATION OP LEADING CHARACTERISTIC LINE 


Input Data 

Calculations of a characteristic network including the input quantities 
when computed by hand calculators are both lengthy and time consuming, particu- 
larly when large numbers of input points along the leading characteristic line 
are desired. A program has therefore been set up to machine compute and deliver 
the input data in punchcard form, which may be placed directly in the charac- 
teristic network program. Through additional processes the data may be indi- 
cated numerically on the cards or tabulated for visual inspection. The fol- 
lowing paragraphs describe the important features of the leading characteristic 
line and give an example of the required input information. 

Regardless of the type of nozzle, either contoured for parallel flow at 
the exit or conical with radial-exit flow, several initial conditions depending 
only on nozzle design are immediately known, namely: Mj, 0^, p., Vjj, and 7. 

With these nozzle characteristics it is possible to determine the leading char- 
acteristic, or Mach, line which separates the nozzle flow from the external 
influence of ambient pressure. The only limitation is that the nozzle must be 
operating in an underexpanded condition such that the pressure ratio P j/Pcc - s 

greater than 1.0. 

Conical nozzle .- The inherent radial flow properties of a conical nozzle 
result in a curved leading characteristic line; therefore, all flow conditions 
along this line are continually varying. It is therefore necessary to deter- 
mine for each x and y location along this line the flow properties required 
in the calculations. It should be pointed out that all phases of the character- 
istic program calculations are performed below the nozzle center line , and that 
all procedures are set up on this basis. 

Point spacings along the leading characteristic are determined, for con- 
venience, in terms of Mach number increments. At low pressure ratios large 
point spacings may be used whereas for higher pressure ratios a much closer 
spacing is required immediately downstream of the nozzle exit. No criterion 
for explicit point spacing, or Mach number increments, exists and, therefore, 
individual judgment must be used depending on the length of the leading char- 
acteristic and the Mach number change from the nozzle exit to the intersection 
of the leading characteristic with the nozzle center line. Several groups of 
Mach number increments are usually calculated; for example , the leading char- 
acteristic shown in sketch (Al) has been divided into three different groups 
of Mach number increments - small increments in group 1 and larger increments 
in groups 2 and 3- If the characteristic program terminates before reaching 
the maximum jet boundary, it may be necessary to use smaller Mach number incre- 
ments in group 1. 
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Center line 



Sketch (Al) 


As an aid in selecting the Mach number increments along the leading char- 
acteristic line, the Mach number at the intersection of the leading character- 
istic with the nozzle axial center line is determined by the method of refer- 
ence 18 . The Prandtl-Meyer expansion angle at the center line is: 


= v N + 29 n 


(Al) 


The change in flow direction between the nozzle exit and end of the leading 
characteristic reduces to 0jg, since the flow direction is zero at the center 
line. The Mach number corresponding to the Prandtl-Meyer expansion angle 

may be obtained from gas flow tables (refs. 29 and 30, for example) at the 
appropriate value of 7. Other parameters needed, in addition to the selected 
Mach number increments, are: Mj , 0 ^, and 7. From the input data the com- 

puter program determines the flow properties at various points along the 
leading characteristic line until the nozzle axial center line is reached. 

Input data downstream of the end of the leading characteristic are not required 
since the characteristic network program (appendix C) contains provisions for 
computing these flow properties. 

The end point of the leading characteristic (on the nozzle center line) 
may be obtained from the following expressions given in reference 18 for cal- 
culating the shape of the leading characteristic. (See sketch (Al) . ) 


R = C 


1 + 


7 



7+1 

M7-1) 


(A 2 ) 


x = R cos 0 - R]\j cos 0 N (A 3 ) 

From the initial flow properties just at the nozzle exit, Mj , , 7, and 0jj, 

the constant C may be determined. By using this constant and the Mach number 
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obtained from expression (Al) and gas flow tables, with its corresponding value 
of (i, the expressions (A2) and (A3) yield the desired x location. 

Other basic equations used in this program, in addition to (A2) and (A3), 
are as follows: 


4 


tan -1 ~*~ 

||m 2 - 1 


(A4) 


V = 



(A5) 


y = 


-R sin 


e 


(A6) 


Example of input data for a conical nozzle.- For a nozzle with the known 
exit conditions of Mj = 5*0, 0]\j = 15°, and 7 = 1.4, the additional required 

engineering input is obtained in the following manner. From the known exit 
conditions and equation (Al) the values of and M at the point of inter- 

section of the leading characteristic with the axial center line are found to 
be = IO6.92 0 and M = 12.02. The Mach number increments are arbitrarily 

chosen as follows: 



^limit 

0.05 

7.0 

.10 

10.0 

.15 

12.1 


The final Mach number of 12.1 was chosen slightly larger than the computed value 
of 12.02 in order to allow the program to determine a more accurate final Mach 
number. Since three different Mach number increments were selected, there will 
be three input cards - each containing the known exit conditions and one of the 
Mach number increments with its corresponding final Mach number. The computer 
program uses the data contained on the three input cards and calculates the 
coordinates and other flow conditions required for the characteristic network 
calculations . 

Contoured nozzle .- In contrast to the conical nozzles, the leading char- 
acteristic associated with contoured nozzles designed for parallel exit flow is 
a straight line as shown in sketch (A2). As a result, all flow conditions 
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along this line are constant. The inclination of this line with respect to 
the nozzle axial center line is equivalent to the Mach angle 4 corresponding 
to the nozzle-exit Mach number. Closely spaced points near the nozzle exit. 


Center line 



Sketch (A2) 

similar to the point spacing discussed for the conical nozzle, are again 
recommended. 

Sonic nozzle .- The leading-characteristic input conditions for a sonic 
nozzle are similar to those for the contoured nozzle. However, if Mj is 

exactly equal to 1.0, the slope of the leading characteristic is infinite and 
the characteristic program will not function under these conditions. For this 
reason, the initial conditions of Mj = 1.0038 and 4 = 85° were selected for 

the current program. Reference 18 showed that a similar selection resulted in 
negligible effects upon the free- jet boundary shape. 



Computer Program 

A computer program has been developed for calculating the leading char- 
acteristic line of a conical nozzle. This program, written in FORTRAN language 
(see ref. 3l)> may be used on the IBM 709^ electronic data processing system. 

The input requirements for this program consist of: (l) the nozzle-exit 

Mach number Mj, (2) the flow direction at the nozzle lip 6^, (3) the ratio of 

specific heats 7, (k) the Mach number increments and (5) the limiting 

Mach number for each Mach number increment for which the flow conditions are 
to be determined. From these input data a series of points defining the 
leading characteristic line are calculated, and cards are punched to be used 
as input for a computer program that calculates the characteristic network of 
a free- jet boundary (see appendix C). 

This particular program is capable of computing the leading characteristic 
for several different Mach number increments during one machine run. The pro- 
gram output for each computed point is x, y, 0, and M and is as follows: 
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P—5432 


leading characteristic LINE program 

PUNCHES INPUT CARDS FOR FREE JET INVESTIGATION (P-5430) 


C 
C 
C 

C INPUT FOR A CASE 

C 

C THETA* MACH NO.. LIMITING MACH NO.* DELTA MACH* GAMMA 

C 

OEG-57 • 295T8 
WRITE (6*50) 

50 FORMAT ( 1 HI /28H0LEADING CHARACTERISTIC LINE » 4X6HP— 5432 ) 

C 

C INPUT A CARD 

C 

7 READ (5*15) FTH .FMACH* EMACH « DEL ♦ GAMMA 
15 FORMAT (5F7.0) 

WRITE ( 6*55 ) GAMMA 
55 FORMAT (///7H0GAMMA* 1 F7 * 3/// 

15X4HMACH16X1HX19X1HY19XI 1HTHETA* 0EG*9X8HMU* DEG*///) 

DU c O .0 
XC*=0*0 
FTH*FTH/DEG 
Cl-GAMMA-1 *0 
C2=GAMMA+1 .0 
C3=FMACH**2 
C4=C3-1 *0 
C5*SQRT(C4 ) 

C EQUATION A4 

FMUC=ATAN( 1 *0/C5 ) 

A-2.0/C1 
B*C2/(4*0*C1 ) 

D* I .0/Cl 
SA«SQRT(C2/C! ) 

SA 1=1 *0/SA 

SB = SA*ATAN(SA1 *C5 > 

C EQUATION A5 

U-SB— AT AN ( C5 ) 

RN«1 • 0/S I N (FTH ) 

C«RN/( ( ( 1 * 0+A* ( 1 . 0/C3 ) )**8 )* (FMACH**0 ) ) 

R*RN 
GO TO 8 

2 I F ( TNET C ) 6 » 7 * 6 
6 FMACH=FMACH+DEL 

C 

C EXIT IF MACH GREATER THAN LIMITING VALUE 

C 

IF (EMACH-FMACH) 7*75*75 
75 C3=FMACH**2 
C4=C3“ 1 .0 
C5*S0RT (C4 ) 

FMUC*ATAN(l *0/C5) 

SB=SA*( ATAN (SA1*C5) ) 

U2-SB— ATAN ( C5 ) 

DU=U2— U 

C EQUATION A 2 

R*C* ( U . 0+A* ( 1 * 0/C3 ) ) *»B ) * (F MACH**D ) 

8 DTH=-.5*DU 
TH2«FTH+DTH 
THETC*-TH2*DEG 
IF (THETC ) 4*7.4 

4 Xl=XC 

C EQUAT ION A I 

U2L*2.0*PTH+U 
C EQUATION A3 

XC*R*C0S(TH2)-(RN*C0S(FTH) ) 

C EQUATION A6 

YC=-R*SIN(TH2 ) 

IF (THETC ) 1 1.11*12 

1 1 U21 *U2 

fmuc=fM^*deg 

PUNCH ANO PRINT 

PUNCH 17*XC.YC*THETC*FMUC 
17 F ORMAT ( 4E 1 6 • 8 ) 

WRITE <6*60 )FMACH«XC ,YC. THETC *FMUC 
60 FORMAT (5E20. 8) 

GO TO 2 

12 UR“<U2L-U21 )/(U2— U21) 

FMACH*UR*DEL+ (FMACH— DEL ) 

XC=UR*(XC-X1 )+Xl 

FMUC=ATAN(1 .0/SQRT(FMACH**2-l .0 ) ) 

THETC*0 • 0 
YC*0.0 
GO TO 1 1 
END 
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APPENDIX B 


DETERMINATION OF CORNER EXPANSION FAN 


Input Data 

An additional input requirement for computing a characteristic network is 
the corner expansion fan just at the nozzle exit. A computer program has been 
developed for a high-speed digital computer to provide these data in suitable 
form to be used directly in the characteristic program. Primary consideration 
was directed toward expediency by limiting the amount of engineering input 
required. 

Two-dimensional flow can be assumed to exist at the nozzle lip. A further 
assumption of instantaneous expansion permits evaluation of all flow conditions 
at the coordinates x / r j = 0 and y/ r j = -1.0. The initial consideration con- 
sists of determining the change in Mach angle \i between that for the design 

exit Mach number and the final value as fixed by the particular pressure ratio 

Poo/P t j P° r w kl c k a characteristic net is to be computed. This change in Mach 

angle ji is then divided into a number of increments, which determines the 

value of Af x to be used in the input data. 

The following table illustrates the general procedure for supplying the 
input data: 


Initial Final 


M Ik 9987008 7.2756624 

e, deg 15.000000 30 . 462664 

deg 11.539999 7.8999996 

v, deg 76.908041 92.370705 

7, deg 1-5999999 


For the known initial conditions at the nozzle exit it is desired to find the 
final conditions for Pj^/Pco = 10-095 in order to compute the input data. From 

the nozzle-exit pressure ratio Pj^Pt j at 7 = 1.4 and the static-pressure 

ratio Pj/P,*, the boundary Mach number corresponding to the pressure ratio 

Poo/Pt , j is found to be 7*2756624. The corresponding values of p and v at 

the boundary Mach number are obtained from flow tables (ref. 28, for example). 
The final flow direction 9final is synonymous with the initial turning angle 

ccjj and is found from the expression 
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0 final = a** = V;L - v N + 0 W (Bl) 

The value of Zp., although selected somewhat arbitrarily, may then be deter- 
mined by the number of rays desired for the variation in 9 . Thus , if 9 
changes by approximately 15-5° and- one expansion ray is desired for each 1 ° 
change in 9, then the change in Mach angle p between the initial and final 
values should be divided into increments of 0.23° (i.e., Zp = 0.23°) • The 
number of significant figures for which Zp is given is unimportant since the 
computer program assumes this is an absolute number and gives results for other 
parameters out to eight significant figures. However, the initial and final 
values of p should be determined as accurately as possible. The required 
input to the computer program utilizing the change in Mach angle p between 
the initial and final conditions is listed as follows: 


x/r^ . . . . 
y/r i .... 
deg . . . 

^initial* de S 
/V y deg • • • 

Hfinal> de S * 


. . 1.4 

. . 0.0 

. . - 1.0 
. . 15.0 

11.539999 

. . 0.23 
7.8999996 


Computer Program 

A computer program has been developed for calculating the comer expansion 
rays to be used as input data in the characteristic network calculation proce- 
dure and is applicable to either conical or contoured nozzles. Input require- 
ments for this program consist of initial flow conditions at the nozzle exit, 
p, 9 jj, and 7 , a value for Mach angle increments Zp, and the value of p 
corresponding to the boundary conditions. The Mach angle increment is deter- 
mined from the fineness desired in the characteristic network. The program 
computes the various flow conditions corresponding to the designated Mach angles 
and cards are punched, which may be placed directly in the characteristic pro- 
gram (appendix C). This program is written in FORTRAN language for the IBM 7094 
electronic data processing system and contains equations (or versions of) (A4) , 
(A5), and (Bl). It is as follows: 
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C P-5433 CORNER EXPANSION RAY PROGRAM 

C PUNCHES INPUT CARDS FOR FREE JET INVESTIGATION ( P-5430) 

C 

C INPUT FOR A CASE 

C 

C 1ST DATA CARO 

C FIRST MU. FIRST THETA. GAMMA 

C DATA CARDS 2.3.4....*. 

C MU VALUE. DELTA MU. STOPPING MU VALUE 

C 

CON- • 1 T453292E— 0 1 
XA-0.0 
YA--1 .0 
WRITE (6.*) 

2 FORMAT! 1 HI /21HQC0RNER EXPANSION RAY4X6HP-5433 ) 

C 

C READ 1ST DATA CARD 

C 

read (5.4 )FMUA .THET A. GAMMA 
4 FORMAT! IF9.0.2F7.0 ) 

THETA--THET A 
WRITE ! 6 . 6 ) GAMMA 
6 FORMAT (///7H0GAMMA-F7. 3/// 

15X4HMACH1 6X1HX19XIHY19X1 1HTHETA. DEG.9X 
28HMU. DEG . I 2X7HV . DEG.///) 

TH2 -THETA *C ON 
FMU2«FMUA*CON 
H- (GAMMA- 1 .0 )/ ! GAMMA* 1 .0 ) 

HI -SORT ! (GAMMA + 1 .0)/ (GAMMA-1 .0 > ) 

C 

C COMPUTE STARTING VALUES 

C EQUATION A4 

C 

FMACH-1 . 0/S I N (FMU2 ) 

ABCD-FM ACM**2— 1 .0 

VPAR-H I *ATAN ! SORT !H*ABCD ) ) 

C EQUATION A5 

V2-VPAR— ATAN !SQRT !ABCD) > 

VPR-V2/CON 

C 

C PUNCH AND PRINT 

C 

PUNCH 10.XA.YA.THETA.FMUA 
10 FORMAT I4E16.8) 

WRITE I6.15JFMACH.XA.YA.THETA.FMUA.VPR 
15 FORMAT ! 6E20. 8 J 
C 

J SAVE THETA AND V 

c 

VI -V2 
TH1 -TH2 
C 

C INPUT DATA CARD 

C 

20 READ ! 5 • 4 ) FMUA . DEL » ENDMU 
30 FMU2»FMUA#CON 

FMACH-1 . 0/S I N !FMU2 ) 

ABCD-FMACH**2-1 *0 

VPAR-H 1 *ATAN ! SORT !H*ABCO ) J 

V2-VPAR— AT AN (SORT ! ABCD) ) 

C EQUATION Bt 

TH2-V1-V2+TH1 
THETA -TH2/C ON 
VPR-V2/C0N 
C 

C PUNCH AND PRINT 

C 

PUNCH 10.XA.YA. THETA, FMUA 

WRITE (6. 15 JFMACH. XA.Y A. THETA. FMUA, VPR 

C 

C SAVE THETA AND V 

C 

VI -V2 
TH1 -TH2 
C 

C DECREMENT MU 

C 

FMUA -FMUA— DEL 
C 

C CHECK STOPPING MU 

C 

I F ! ENDMU— FMUA J 30 , 30 . 20 
END 



APPENDIX C 


CHARACTERISTIC METHOD FOR CALCULATING FREE -JET BOUNDARIES 


General Description of Characteristic Network 

The calculation of free -jet boundaries by the method of characteristics 
makes use of three-dimensional irrotational equations of flow. The character- 
istic net describing the flow field, a lattice point-type structure requiring 
no iteration, is machine computed in a point-to-point calculation procedure. 
Computations begin with a set of given conditions on the leading characteristic 
plus the conditions for a series of two-dimensional expansion rays originating 
at the nozzle lip. The existing program, permitting calculations to high pres- 
sure ratios and large distances downstream of the nozzle exit, is the result of 
extensive modifications and refinement of the program used in reference l 8 . 

The general assumptions concerning the characteristic network are: 

(1) The characteristic solution assumes a constant ambient pressure and 
hence constant Mach number boundary. 

(2) The ratio of specific heats 7 remains constant in all calculations. 

(3) Isentropic flow is assumed to exist throughout the flow field. 

Various sections of a typical characteristic net require different equa- 
tions to determine the local flow properties. These sections are listed as 
follows in order of their probable occurrence during the calculation procedure: 

(1) Solutions for general points. 

(2) Solutions for points on the jet boundary. 

(3) Solution for the intersection of characteristic lines of the same 
family indicating the presence of an internal shock. 

(k) Solution for a point adjacent to the nozzle axial center line. 

(5) Solution for points on the center line. 

Each of these topics is discussed in the following sections. 
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Sketch (Cl) 

The types of points in sketch (Cl) are identified as follows: 

(1) nozzle lip, source for a series of given two-dimensional expansion rays 

( 2 ) ,(3),( i l-),(5) input data points on leading characteristic line 

(6) last leading characteristic input data pofht (also first point on 

center line) 

(7) >(8),(9),(l0),(ll),(l2),(l3) general points 

(14) boundary points 

( 15 ) same family points (internal shock) 
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( 16 ) points adjacent to center line 

( 17 ) center -line point 


Symbols 

M Mach number 

T temperature, °R 

V velocity, ft/ sec 

Vz limiting velocity, ft/sec 

W ratio of local velocity to limiting velocity 

x, y Cartesian coordinates, referred to nozzle radius 

7 ratio of specific heats 

0 flow angle between velocity vector and X-axis (positive in first quad- 

rant and negative in fourth quadrant), radians 

p Mach angle, sin"l i, radians 

M 

p density, slugs/cu ft 

Subscripts: 

A points of first family 

B points of second family 

C points calculated from A and B 

t stagnation conditions 


General-Point Solution 

In supersonic flow if the flow conditions at two points near each other are 
known, the flow conditions at 4 some point C, which is the intersection of two 
characteristic lines of different families starting from the two known points, 
can be determined. In sketch (Cl) calculations begin with the given conditions 
at points (l) and (2) and proceed along the line containing points (2) and (l4) 
with use of the general-point solutions to find points ( 7 ) to (11) and the 
boundary -point solution to find point (l4) . The procedure is then repeated 
beginning with the given point ( 3 ) and the conditions as computed for point ( 7 ) 
to find point (12) ; the points (12) and (8) are used to find ( 13 ) , and so forth. 
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An enlargement of points ( 3 ) and ( 7 ) as used in obtaining the general solu- 
tion for point (12) is shown in sketch (C2). 


X A X B X C 



With all conditions known for points A and B, the flow conditions at point C 
remain to be found. The upward sloping line AC, a characteristic of the first 
family, is defined as 


dy _ y C ~ y A 
dx Xq - x^ 


tanfe. + 



(Cl) 


The slope of line BC , a second-family characteristic, is defined as 


dy _ y C ~ y B 
dx xq - x-g 


tan^ 0 B - u B ) 


A simultaneous solution of equations (Cl) and (C2) yields 


Xc = 


x A tan^e A + 

4 A ) 

1 

+ 

£ 

- x-g tan | 

( 0 B “ 4 B ) 

tan| 

( e A 

+ 4a) 

1 - tan j 

( 0 B “ 4b) 

1 


(C2) 


(C3) 
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Now that xq is known, equation (Cl) may he solved for the y coordinate, 
y c , from 


yc = (*c - x A ) tan ( e A + ^a) + y A (CM 

First- and second-family equations in terms of velocity are as follows (from 
ref. 32, p. 2 6 k): 

For the first family 


&a ^A 

— d0 * tan (j. . — Z« — — — 0 

V A A A y A 


For the second family 


dV B 

— + d0 B tan 


dx-o 

'V'" 


(C5) 


(C 6) 


For definitions of Z A and mg, see equations (C9) and (CIO). The velocity 

terms may he reduced to the nondimens ional ratio form of W hy dividing through 
hy the limiting velocity Vj. 



Therefore, in equations (C5) and (c6) 



W C - W A 



w c - w B 


and 


de B = de A + e A - 0 B 


Substituting these expressions into equations (C5) and (C 6) and solving simulta- 
neously results in 
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|tan [i-g + 


W, \ W« / \ dx B dx A W A 

^ - 1 - ^ - ( 9a ' OBjtan , B + ^ J - Z A — gr (C7) 


which is similar to the form in reference 32. By substituting 


d0A = (Qq - 0a) 


and 


= (xq “ x a) 


= (*C “ X B) 


equation (C7) is reduced to the following form with only one unknown: 

-0 A tan h a + “ x a)J + W B + w b|®B tan ^B + ^(*0 " x b)| 


-% - w A 


0 C = 


(w A tan HA + W B tan Hb) 


where 


*A = 


sin h a sin 0 A tan h a 


cos 


( 0 A + ^a) 


and 


m B = 


sin h b sin 0 B tan h B 


COS I 


( e B “ ds) 

Now that 0 q is known, equation (C5) is used to obtain 


W C = w A + W A 


tan h a (0 c - 0 A ) + ^(xc - x a) 




The remaining parameters to be determined at point C are 


^ 


(c8) 


(C9) 


(CIO) 


(Cll) 


(C12) 
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1 


Me = 


sin 




(C13) 


and 



(Clk) 


(C15) 


Jet Boundary Point 

Boundary-point calculations, although still dependent upon the basic first- 
and second- family relationship, require a different form of the equations for 
determining the conditions and location of the unknown boundary point C. As in 
previous calculations it is necessary to adhere strictly to the signs associated 
with the various angles. For these calculations (see sketch (C3)) it is noted 



Sketch (C3) 


30 


that the line AC is actually a streamline, where previously the line AC has 
been a characteristic line. The slopes of the streamline AC (first family) and 
characteristic line BC (second family) are 


dy o y C - y A = tan (cl6) 

dx Xq - x A 


dy _ y C ' y B 
dx Xq - x B 


tan 


(®B " ^b) 


(C17) 


The x location as determined from a simultaneous solution of equations (Cl6) 
and (C17) is 


Xn = 


x A tan e A - y A + y B 


x E tan(8 E 


~ ^b) 


tan 0 A - tan^0 B - |i B ^ 


(cifl) 


and from equation (Cl6) 


y C = (*C 


x a) tan 0 A + y A 


(C19) 


Since the boundary Mach number is constant, 


[i c = n A and W c = W A ; the flow 


direction at point C is obtained from the basic second-family form of equa- 
tion (C6). In equation (c6) 


de B = 0c - ©B 


and 


W C ■ W B *C - X B 

+ m B 


e c = 0 b + 


Wb 


yb 




tan [ x-q 

where m-g is as defined in equation (CIO) . 


(C20) 


Intersection of Same Family Characteristic Lines 

The presence of a finite ambient pressure requires a turning of the jet 
boundary in a downstream direction which in turn introduces a series of weak 
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compression waves. These waves,, which are shown in sketch (c4) and are members 
of the first-family type of characteristic lines (see sketch (Cl)), eventually 
intersect with the expansion rays originating at the nozzle exit or with each 
other and form the internal shock. The intersecting characteristic lines, both 
members of the first-family type of lines, are shown in sketch (c4). Prior to 
the occurrence of this intersection at point (15), calculations progressing 
along line 7 in sketch (c4) (arbitrarily identified lines and points) have 
involved members of different families. The normal procedure of using points 
(4) and (13) as A and B, respectively, to find (l4) as C, and then using ( 5 ) 
and (l4) as new values of A and B to find point (l6), results in a foldback 


^o 



Sketch (C4) 

of the characteristic net (see ref. 18). Once the foldback has begun it becomes 
progressively more severe and eventually a point is reached where the program 
can no longer continue. In order to eliminate the foldback, flow conditions 
must be found for the point of intersection ( 15) ; once these are found, then 
points (15) and (13) may be used as A and B points, respectively, in the 
general-point equations . 

An enlargement of the two intersecting rays is shown in sketch (C5) • Since 
(i is always considered to be positive and strict observance must be made of 
the sign associated with the angle 9, the characteristic lines are both members 
of the first family and the slopes are 

/dy\ = yc - ^A 
WA *C - X A 


tan(9 + p); therefore, 

= tan(9 A + pa) (C21) 
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Simultaneous solution of these two equations yields 



tan 


( e A + Ha) - y A + y B - x b tan ( e B + » b ) 
tan^0 A + - taring + ^b) 


(C23) 


and from equation (C2l) y q is found to he identical to equation (C4). 

In obtaining 9 q at the intersection of two rays of the same family, it 

is to be remembered that the primary difference between the general-point and 
same-family solutions is a change in notation of the angle which the character- 
istic lines make with respect to the X-axis - (0 + p) for the first family and 
(0 - p) for the second family. Since the intersection of rays of the first 
family is involved, (0 + p) must appear in all equations in which there is a com- 
bination of 0 and p. It is therefore possible to use the general-point solu- 
tion for 0 q with the exception that the m^ term must include cos(0 + p). 
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The equation for Qq is the same as equation (C8) except that 


m B 




sin |a-g sin 8-g 
cos^Gb + Hb) 


(C2k) 


All other terms are identical to those of the general -point solution. 


Point Adjacent to the Center Line 

The characteristic lines used in defining flow conditions of a point adja- 
cent to the axial center line involve a point on the center line, point B in 
sketch (C 6), for which the flow direction is 0g = 0°. The general -point equa- 
tions (C3) and (c4) provide solutions for hoth x and y coordinates of 


B Vg Nozzle center line 



point C; however, complications arise in determining 0q since the term 
cLx 

mg — used in the second-family equation (C 6) becomes indeterminate at the 


axis. Reference 32 makes several simplifying assumptions and derives an equa- 
tion for the second- family characteristic line BC which is: 


W B 


(C25) • 
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Equation (C25), written for a point above the center line, conforms to condi- 
tions below the center line in this report by changing the minus (-) sign in 
front of the term ^A0 tan to a plus (+) sign. An equation now exists for 

a second-family characteristic line similar in form to the second-f ami ly line 
in equation (C 6). However, in equation (C25) 


and 


Ax _ X C - *B _ 1 

yc tan |i B 


sin p B 
cos 


tan 


Substitution of these terms into equation (C25) along with the required sign 
change produces 


W C - W B 
W B 



0 B )tan n B = 0 


(C26) 


The simultaneous solution of the first-family equation ( C5) and the second- 
family equation (C2 6 ) yields 


" W A “ w a(- 0 A tan 4A + l A X °y A —j + W B + sw B e B tan Hb 
6 C = A “ TT *: -f- 


W A tan + 2W B tan p B 


(C27) 


The velocity ratio Wq and Mach angle |J.q may now he found by using the 
general -point solutions. 


Center -Line Point 

The procedure for finding flow conditions at points along the nozzle axial 
center line is simplified in that both y^ and 0 q are zero. The intersection 

of the leading characteristic with the center line is the last given point in 
the input data; therefore, the intersection of the first expansion ray with the 
center line is the first point of this type which must be calculated and is 
shown in sketch (C7). An identical expansion originating at the opposite side 
of the nozzle results in point B which is a mirror image of point A. As a 
result of this similarity only the first-family characteristic line AC is needed 
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Sketch (CT) 

to determine the location of point C. With y^ and 0 q equal to zero equa- 
tion (Cl) reduces to 

y A 

xc = x A (C28) 

tan^0 A + |i A j 

Other parameters such as Wq, hq, and Mq may he found with the general -point 
equations. 


Computer Program 

A computer program has been developed for calculating the characteristic 
network of a jet boundary using three-dimensional irrotational equations of 
flow. This program is written in FORTRAN language (ref. 31) and is presented 
at the end of this appendix. 

A set of initial conditions, x, y, 0, and \i, defining a series of 
points along the leading characteristic line is supplied to the computer in 
addition to the conditions defining the point origins of a series of two- 
dimensional expansion rays originating at the nozzle lip. The number of rays 
given, for any given exit Mach number, depends upon the fineness of the charac- 
teristic net desired and the pressure ratio in question. Computer programs have 
been written that may be used to supply the initial conditions for the leading 
characteristic line (see appendix A) and the conditions for the corner rays 
( see appendix B) . 


36 



Points are computed along a downward sloping characteristic line from each 
given point on the leading characteristic line until the boundary point is com- 
puted. When the leading characteristic line is ended, points are computed 
along the nozzle axial center line with two new points being computed by inter- 
polation between normal successive center-line points. Additional character- 
istic lines are computed originating at these center-line points. When |y| 
has reached a maximum on the boundary, computation of points along character- 
istic lines continues until x on the center line or on the leading character- 
istic line has reached the value of x at the maximum boundary. 

The following program has been used on the IBM 709^ electronic data proc- 
essing system at the Langley Research Center to obtain the results presented 
herein. 
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P—5430 FREE JET INVESTIGATION - CHARACTERISTIC METHOO 

INPUT FOR A CASE 

1ST CARD IDENTIFICATION* COL. 1-8© 

2ND CARD NO. OF CORNER RAY CARDS. COL. I -3 

NO. OF LEAOING CHARACTERISTIC CARDS* COL.4-6 
VALUE OF LAST Y (RAY ) READ IN. COL. 7-1 3 
GAMMA. COL. 14-20 
DATA CAROS 3*4. ....... 

1ST GROUP 

EACH CARD (NO. OF CAROS GIVEN BY COL.1-3 ON CARO 2) 
DESCRIBES A POINT ON THE CORNER RAY ANO CONTAINS 
4 VALUES ( X *-Y* -THETA, MU > 

2ND GROUP 

EACH CARD (NO. OF CAROS GIVEN BY COL. 4-6 ON CARD 2) 
DESCRIBES A POINT ON THE LEADING CHARACTERISTIC 
LINE AND CONTAINS 4 VALUES (X.-Y. -THETA* MU) 

DIMENSION CASE ( I 4 ) » XA ( 400 ) . YA (400 > .TMETA (400 ) .FMUA (400 ) « 
LWA (400) .CHARAC (6.700 ) .XC (400 ) .YC (400 ) .THETC UOO) . 
2FMUC(4<10) .WC<400> .FMACHC ( 400 > 

COMMON CASE. XA.YA. THETA. FMUA.WA, 

1 FMACHC . XC . YC . THETC . FMUC . WC . CHARAC . 

2GM1 .GMI2.XP.LINE.NA.NC.NRAY.NREFL 
C 

C INITIALIZATION SECTION 

C 

c 

C EQUATION (Cl 3 ) 

FMMU(X>«1 ./-SIN(X) 

DEG-57. 29578 
WRITE (6.1) 

1 FORMAT! IH1///23H0FREC JET INVEST I GAT I ON4X6HP-5430//V ) 

2 CALL SLITE (0) 

nrefl*o 

CALL DVCHK (KOOOFX) 

GO TO( 4 .4 ) .KOOOFX 

c 

c READ IN CASE 

C READ 1ST 2 CAROS 

C 

4 READ (5*5 )CASE *L INE.NCHAR* YPREV .GAMMA 

5 FORMAT ( I3A6.A2/2I3.2F7.0) 

NRAY-LINE-1 

C 

Q READ CORNER RAY CARDS 

C 

READ (5.7) (XC ( I ) . YC ( I ) .THETC t I ) .FMUC ( I) . I *1 .LINE) 

7 FORMAT (4E1 6.0) 

C 

C READ LEADING CHARACTERISTIC CAROS 

C 

READ (5.7) ( (CHARAC! I. J).I»1 .4 ) . J» 1 , NCMAR ) 

C 

C PRINT ID AND GAMMA 

C 

WRITE (6. 9} CASE .GAMMA 
9 FORMAT (1XI3A6«A2^/7M0GAMMA»F7.3///) 

C 

C SET READ TO DETERMINE LAST CHARACTERISTIC POINT GIVEN 

C 

I-NCHAR-I 
DO 10 J-t.I 

10 CHARAC(6. Jl-0.0 
CHARAC ( 6. NCHAR>»1 .0 
WR I TE (6.11 > 

11 F0(RMAT(6H0 MACH1 3X1 HX1 6X1 HY1 6XSNTHETA 

I 1 1 X2HMU I 4X1 HWt SX4HT/T0 1 2X8HRHO/RHOO// ) 

c 

GM1 -GAMMA- 1 .0 
GMI2-GM1/2.0 
XP-1 .0/GM1 
C 

C CHANGE THETA ANO MU TO RADIANS 

C COMPUTE WC ANO MACHC FROM INPUT MUC 

C 

DO 19 I - 1 .LINE 
THETC ( I ) -THETC ( I )/OEG 
FMUC 1 1 ) -FMUC ( 1 1 /DEG 
C EQUATION (Cl 2 ) 

WC(1 )-SORT(GMl/(2.*SIN(FMUCtI))**24GMI )) 

FMACHC ( I ) -FMMU (FMUC ( I ) ) 

19 CONTINUE 
C 

C PRINT CORNER RAY 

C 

DO 22 I-I .LINE 

CALL PRINT (FMACHC (I ).XC( 1 ) *YC( I ) * THETC ( I 1. 

1FMUC ( I I.WCttn 
22 CONTINUE 
CALL PREND 
C 

C MOVE C ARRAY 

C 

DO 23 I - ! .LINE 
XA ( I ) -XC ( I ) 

y«ut*.vc4i> 



THeTA<t»*?Hfcrem 

FMUA ( I ) -FMUC ( l ) 

23 wAm-wcm 

c 

C CHANGE CHARACTERISTIC THETA AND MU TO RADIANS 

C COMPUTE ALL W FOR LEADING CHARACTERISTICS 

C 

DO 26 1-1 tNCHAR 

24 CHARACI3* I )*CHARAC(3*I >/DEG 
CHAR AC (4.1) -CHARAC ( 4 • ! )/DEG 

C EQUATION (Cl 2 > 

CHARAOS.I )«SQRT(GM!/(2**SIN(CHARAC(4. m**2+«Ml )> 
26 CONTINUE 
C 

C COMPUTE 1ST LINE-NO REFLECTION 

C 

28 I CELL-2 
LEAOC-1 
LINE-NRAY-1 
NA-1 
NC- I 
C 

C MOVE 1ST LEADING CHARACTERISTIC TO 1ST C POINT 

C 

32 XC-CHARAC (1*1) 

YC-CHARAC (2.1) 

TH£TC«CHARAC(3. 1 > 

FMUC-CHARAC(4. 1 ) 

WC-CMARAC(S«1 > 

WRITE (6* II) 

FMACHC -FMMU ( FMUC 1 

CALL PR TNT ( FMA CMC . XC ♦ YC . THET C * FMUC . WC ) 

35 NA-NA+t 
NC-NC+l 

CALL GENL (XA ( NA > . YA ( NA > , THET A ( NA ) .FMUA ( NA ) « WA < NA ) . 

1 XC (NC— 1 ) » YC(NC— I ) *THETC(NC—1 >»FMUC(NC~1 )«WC(NC-1 >. 
2XC (NCI. YC(NC). THETC (NC) .FMUC (NC).WC(NC) ) 
FMACHC(NC)«FMMU(FMUC (NCI ) 

CALL PR I NT (FMACHC (NC 1 « XC (NC ) * YC (NC ) tTHETC (NC > . 
1FMUC(NC).WC(NC) ) 

LINEaLINE-I 
IP (LINE >998.40.35 
40 NA-NA+1 
NC-NC+1 

CALL BNO«Y(XA(NA) «YA (NA) .THETA (NA )* FMUA ( NA ) .WA(NA) • 
1 XC ( NC— 1 ) * YC (NC— 1 >.THETC(NC-1 ).FMUC(NC-1 ) .WC(NC-1 )« 
2XC (NC I ♦ YC (NC) * THETC (NC ) «FMUC (NC) ,WC(NC) ) 

FMACHC (NC) -FMMU (FMUC (NC) ) 

CALL PR 1 NT (FMACHC (NC) *XC(NC> *YC (NCI* THETC (NC > . 

1FMUC (NC > * WC (NC ) I 

SIFSSB 1 

C 

C MOVE C ARRAY 

C 

DO 45 I-1.L1NE 
XA (!) -XC ( 1 ) 

YA ( I )— YC ( I ) 

THETA ( I ) -THETC ( I ) 

FMUA (!) -FMUC ( 1 ) 

45 WA ( I )— WC ( I ) 

NREFL-1 

C 

C 

c MODE I - SL2 OFF* READ-0 

C MAX BOUNDARY NOT REACHED* USE CHARAC LINE 

C 

C 

100 MOOE-1 

102 LINE-NRAY+NREFL-1 

leadc-leadc+i 

NA — 1 
NC-1 

C 

C MOVE LEADING CHARACTERISTIC TO 1ST C POINT 

C 

104 XC-CHARAC ( 1 « LEADC ) 

YC-CHARAC ( 2 * LEADC ) 

THETC-CHARAC <3 .LEADC ) 

FMUC-CHARAC (4.LEADC) 

WC-CMARAC (5.LEADC ) 

108 RE AO-CH ARA C ( 6 » LEADC ) 

WR I TE (6.11 ) 

FMACHC -FMMU ( FMUC ) 

CALL PR I NT ( FMACHC . XC . YC * THETC . FMUC . WC 1 

C 

C TEST FOP 1ST CENTERLINE 

C 

I F ( READ ) 998 .110.109 

109 NA—2 
NC-2 

CALL OFCNT (XA (2 ) « YA (2 ) .THETA (2 > « FMUA (2 ) *WA (2 > . 
1XC*FMUC«WC»XC(2 ).YC (2). THETC (2) .FMUC (2) »WC(2> ) 
FMACHC ( 2 > -FMMU ( FMUC ( 2 ) ) 

CALL PRINT (FMACHC (2 ).XC( 2 ).YC(2).THETC(2). 

1 FMUC ( 2 ) . WC ( 2 ) ) 

LINE -LINE-1 

110 NA-NA+I 
NC-NC+t 

CALL GENL CXA (NA ) « YA (NA ) .THETA (NA ) . FMUA ( NA ) * WA ( NA ) » 
1XC<NC-1 ).YC(NC-I I.THETC(NC-1 > .FMUC ( NC-1 ) . WC (NC-1 ). 
2XC«NC).YC (NC). THETC (NO .FMUC <NCI .WCINCIl 



lPtLlNe-4>U5*ii2Uii 

111 CALL test 
NA=NA 

112 CALL TEST2 (ISIG) 

NC-NC 

113 GO TO (114.123). ISIG 

114 FMACHC (NC)-FMMU (FMUC (NC) ) 

CALL PR 1 NT < FMACHC (NC) *XC(NC) . YC ( NC ) * THETC INC)* 

1FMUC (NC) « WC ( NC > ) 

LINE-LINE-1 
IF (LINE) 998. 120* 1 1 0 
120 NC-NC+1 

122 CALL BNDRY(XA (NA) .YA (NA > .THETA (NA ) .FMUA(NA) .WA(NA) » 
1XC (NC— l ) * YC ( NC— 1 ) . THETC (NC-1 ) .FMUC (NC-1 ) .WC(NC-1 ) * 
2XC (NC 1 » YC (NC) . THETC (NC) .FMUC (NC) .WC (NC) ) 

FMACHC (NC)-FMMU (FMUC (NCl) 

CALL PRINT (FMACHC (NC J *XC (NC ) * YC (NC ) * THETC (NC ) * 

1 FMUC (NC > «WC(NC > ) 

LINE-NRAY+NREFL+1 
CALL PREND 

MOVE ARRAYS 

00 125 I-l.LINE 
XA ( I )-XC ( I ) 

YA ( 1 )«YC( I ) 

THETA ( I ) -THETC ( I ) 

FMUA ( I ) -FMUC < I ) 

WA ( I )*WC ( I ) 

CONTINUE 

TEST FOR NEXT MODE 

T F ( YPREV— YC ( NC ) ) 130 *150 .150 
CALL SLITE (2) 

XMAX-XC (NC ) 

IF (READ) 996 .300 . 1 32 

1 CELL-3 
GO TO 400 
YPREV-YC(NC) 

IF (READ )998 * 152*151 
I CELL-3 
GO TO 200 
NREFL-NREFL+1 
GO TO 102 


MODE II - SL2 OFF* READ NOT 0 

MAX BOUNDARY NOT REACHED* CENTERLINE POINTS COMPUTED 


mODE-5 

GO TO (760.740.700) * TCELL 
L I NE-NR AY+NREFL— 2 
NA=NA+1 
NC-NC+1 

CALL GENL (XA (NA ) , YA (NA ) * THET A (NA ) * FMUA ( NA ) * If A ( NA ) * 

I XC (NC— I > . YC ( NC— 1 ) «THETC(NC-1 ) *FMUC (NC-1 ) * WC1NC-1 ) * 
2XC (NC ) . YC (NC) .THETC (NC ) .FMUC (NC ) * WC (NO ) 

IF (LINE-4 >206.206.205 

205 CALL TEST 

na =na 

206 CALL TEST2 (ISIG) 

NC-NC 

207 GO TO (208.212) . ISIG 

208 FMACHC (NC )«FMMU (FMUC (NC) ) 

CALL PRINT (FMACHC (NC ) *XC (NC ) * YC ( NC ) * THETC (NC ) * 

I FMUC (NC ) « WC (NC > ) 

LINE-LINE-1 
IF (LINE )998. 21 0.204 
210 NC-NC+1 

CALL BNORYIXA (NA ) * YA (NA) .THETA (NA) * FMUA (NA ) » WA (NA ) * 
1 XC ( NC— I ) *YC (NC— 1 ) .THETC(NC-1 ). FMUC (NC-1 > *WC(NC-1 ) . 
2XC (NC ) . YC (NC) .THETC (NC ) *FMUC (NC > .WC(NC) ) 

FMACHC (NC)-FMMU (FMUC (NC) ) 

CALL PR I NT (FMACHC (NC ).XC(NC) *YC(NC ) .THETC (NC ) « 

1FMUC (NC ) «WC(NC ) ) 

212 L INE-NRAY+NREFL+1 
CALL PREND 
C 

C MOVE ARRAYS 

C 

DO 215 1-1 .LINE 
XA ( I )«XC( I ) 

YA ( I ) - YC ( I ) 

THETA ( I ) -THETC ( 1 ) 

FMUA (I ) —FMUC ( I ) 

WA ( ! )*WC( I ) 

215 CONTINUE 
C 

C TEST FOR NEXT MODE 

C 

IF(YPREV-YC (NC) >230*220.220 
220 YPREV-YC(NC) 

GO TO 201 

230 CALL SLITE (2) 

XMAX-XC ( NC ) 

GO TO 400 
C 

c 

C MODE III - SL2 ON, READ-0 

C flAX BOUNDARY REACHED. USE CHARAC LINE 


132 

150 


500 

201 

202 

204 



e 

c 

300 MODE -3 

302 Lf NE-NRAY+NREFL-1 

leaoc*leaoc+i 

NA-1 

NC-1 

C 

C MOVE CHAPA C TO 1ST C POINT 

C 

30* XC-CHARACU .LEADC) 

YC-CHARAC { 2 , LEADC ) 

THETC-CHARAC ( 3 .UEAOC ) 

FMUC -CHARAC ( * « LEADC > 

WC-CHARAC (3. LEADC ) 

READ-CHARAC (6.LEADC ) 

IF (XC-XMAX ) 303 , 303.990 

305 WR!TE(6« 1 1 ) 
fmachc -fmmu(fmuc ) 

CALL PRINT ( FMA CMC . XC . YC . THETC » FMUC , WC > 

C 

C TEST FOR 1ST CENTERLINE 

C 

IF (READ 1998. 31 0 * 306 

306 NA-3. 

NC-2 

CALL 0FCNTCXA(2).YAf21.THCTA (21. FMU* (2)»WA(2). 
IXC.FMUC.WC.XCf2 >• YC(2 > .THETC (2 > .FMUC (2 > .WC (2 > > 

IF (XC f 2 ) -XMAX ) 307 .307, 789 

307 FMACHC(2)»FMMU(FMUC(2> > 

CALL PR I NT (FMACHC (21 «XC (2 1 ♦ YC (2 ) .THETC (2) * 

1 FMUC < 2 1 ♦ WC (21) 

LINE-LINE-1 
310 NA-NA+1 
NC-NC+1 

CALL <5ENL(XA (NA >. YA < NA >. THETA (NA ) , FMUA (NA 1 • WA (NA 1 . 

1 XC (NC— 1 1 »YC (NC— 1 >.THETC(NC-I > .FMUC (NC— 1 1 • WC (NC-1 1. 

2XC (NC ) » YC <NC ) . THETC (NC), FMUC < NC > .WC(NC>> 

C 

C TEST FOR MAX BOUNDARY X 

C 

I F ( XC ( NC > -XMAX ) 320 , 320 ♦ 350 
320 TF(LINE-A >330.330,329 

329 CALL TEST 
NA-NA 

330 CALL TEST2 (ISI<5> 

NC-NC 

331 GO TO (333.370). ISIG 
333 FMACHC (NC)«FMMU< FMUC (NC) 1 

CALL PRINT(FMACHC(NC).XC(NC) .YC (NC ), THETC (NC » . 

I FMUC (NC ) ♦ WC (NC 1 > 

LIHe-LI^E-I 
IF (LINE 1998 .370,310 
C 

C COMPUTE RAYS AND REFLECTIONS 

C 

350 LOOK -NC— I 

IF (LOOK-NRAY) 352.333.356 

352 NRA Y=LOOK 

353 NREFL-0 
GO TO 370 

356 NREFL-LOOK-NRAY 
370 LINE-NRAY+NREFL 
CALL PREND 

C 

C EXIT TEST 

C 

IF < NRAY—I >990.990.375 
C 

C MOVE ARRAYS 

C 

375 DO 380 I-l.LINE 
XA ( I ) — XC ( 1 > 

YA ( I )-YC( I ) 

THET A ( I J-THETC ( 1 ) 

FMUA (I 1 -FMUC ( 1 ) 

WA ( 1 )-WC( I) 

380 CONTINUE 
C 

C TEST FOR NEXT MODE 

C 

1 F (READ ) 998 , 302 . 382 
382 I CELL-3 
C 

c 

c MOOE IV - SL2 ON. READ NOT 0 

C MAX BOUNDARY REACHED. CENTERLINE POINTS COMPUTED 

C 

C 

*00 MODE-* 

*02 GO TO (760.7*0.700) , ICELL 
*0* LINE-NRAY+NREFL-2 

IF (LINE) 998, 470. *08 
*08 NA C NA+ 1 
NC-NC+1 

CALL GENL (XA (NA ) . YA (NA ) , THETA < NA > .FMUA (NA > , WA ( NA ) , 

1 XC (NC— 1 1 » YC(NC— 1 ) ,THETC(NC-1 ) .FMUCCNC-l 1 ,WC(NC-1 ) « 

2XC ( NC ) . YC (NC) .THETC (NC), FMUC (NCI .WC(NC) ) 

C 

C TEST FOR MAX BOUNDARY X 

C 

IF (XC (NC l-XMAX )*30.*20.*9« 



*20 1F<L!N£-*}*30.A30.*2$ 

*29 CALL TEST 
NA-NA 

430 CALL TEST2 ClSIG) 

NC-NC 

*3t GO TO (**5.*70>. ISIG 
445 FMACHC (NC)-FMMU (FMUC (NC ) > 

CALL PRINT(FMACMC(NC).XC(NC> ,YC(NC) .THETC (NC>. 

1 FMUC ( NC ) « WC ( NC ) ) 

LINE -LINE-! 

IF (LINE ) 990 *470 • *06 
C 

C COMPUTE pays AND REFLECTIONS 

c 

450 LOOK-NC-1 

IF (LOOK— NR AY )4S2 *453 «4S6 

452 NRAY-LOOK 

453 NREFL-0 
GO TO *70 

456 NREFL-LOOK-NRAY 
470 LINE-NRAY+NREFL 
CALL PREND 
C 

C MOVE ARRAYS 

C 

475 00 450 r*I .LINE 
XA ( I >«XC( I> 

YA(I>«YC(I> 

THETA (! )«TMETC( I) 

FMUA ( I ) -FMUC ( I » 

WAd)-WC(I) 

450 CONTINUE 
GO TO 402 
C 

c 

c interpolation section for centerline points 

c 

C FOR I CELL- 3 

C 

700 I CELL -2 
NA-2 
NC-2 

CALL CENTL (XA (2 ) • YA (2 ) .THETA (2 > .FMUA (2 ) ,WA (2 > . 
1TX.TMU.TW) 

C 

C COMPUTE DELTAS AND SAVE 

C 

DELX*(TX-XA)/3. 

OELMU- < TMU-FMUA ) /3 • 

DELW- (TV-WA >/3, 

C 

C COMPOTE 1^3 DISTANCE FOR Ntw C 

C 

XC-DELX+XA 

FMUC-DELMU+FMUA 

WC-DELW+WA 

C 

C SENSE LIGHT TEST 

C 


712 CALL SLITET(2,K000FX> 

GO TO(7l5.720) *XOOOFX 
715 CALL SLITE (2) 

IF (XC— XMAX ) 720 • 720*990 
720 WR I TE ( 6. 1 1 ) 

FMACHC-FMMU ( FMUC } 

CALL PR I NT ( FMACHC .XC . YC . THET C . FMUC . WC > 

C 

C COMPUTE 1 PT OFF CENTERLINE 

C 

CALL OFCNT (XA (2 ) . YA ( 2 ) « THETA (2 > . FMUA ( 2 > . WA (2 > . 
1 XC . FMUC . WC * 

2XC (2 ) . YC (2 ) • THETC ( 2 ) .FMUC (2 I .WC (2 ) ) 

CALL SLITET (2.KOOOFX ) 

GO TO ( 725.730 ) .KOOOFX 
725 CALL SLITE (2) 

IF (XC (2 I-XMAX >727.727,728 

727 FMACHC (2 ) -FMMU (FMUC (2 > I 

CALL PR ! NT (FMACHC (2>.XC(2).YC<2>. THETC < 2 ) • 

I FMUC (21. WC (2 J 1 
GO TO 404 

728 NR AY- 1 
NREFL-0 
LINE-1 
CALL PREND 
GO TO 990 

730 FMACHC (2 ) mFMMU (FMUC (2 ) ) 

CALL PRINT (FMACHC (2 > ,XC <2 ) .YC(2> .THETC (2), 

1 FMUC ( 2 ) . WC ( 2 ) ) 

NREFL-NREFL+ I 
GO TO 202 

0 

C FOR I CELL -2 

C 

740 I CELL- I 
NA-2 
NC-2 

XC-XA+OELX 
FMUC-FMUA+DELMU 
WC-WA+OELW 
GO TO 712 
C 
6 


FOR I CELL- 1 



e 

760 I CELL-3 
NA-3 
NC-2 

CALL CENTL (XA!2>»YA(2) .THETA 12 J.FMUA (2) .WA (2 ) % 
1XC.FMUC.WC) 

CALL SLITET12.KOOOFX) 

60 T01770.780) .KOOOFX 
770 CALL SL1TE (2) 

IE fXC-XMAX >760.780. 990 
780 WRlTE16.lt ) 

FMACHC -FMMU c fmuc ) 

CALL PRINT1FMACHC.XC.YC.THETC.FMUC.WC> 

CALL OFCNT 1XA 1 3 ) • YA ( 3 ) . THETA ( 3 ) . FMUA 1 3 > * WA 1 3 > • 

1 XC.FMUC . WC.XC 12 ) • YC 12 > • THETC 12 ) .FMUC 12 ) . WC 121) 
CALL SLITET12.KT000FX1 
60 TO 1 785.790 ) .KOOOFX 
785 CALL SL1TE (2) 

T F 1 XC 1 2 ) -XM A X ) 787 . 787 . 789 
787 FMACHC12).FMMU1FMUC12> ) 

CALL PRINT ( FMACHC 12).XC(2) «YC(2) . THETC 1 2 ) . 

1 FMUC 12). WC 12) ) 

NRAY-NRAY-l 
60 TO 404 

789 NRAY-1 
nrefl-o 
LINT-1 
CALL PREND 
60 TO 990 

790 NREFL-NREFL+1 

FMACHC 1 2 ) -FMMU 1FMUC (21) 

CALL PRINT! FMACHC 12) .XC12) .YC12) .THETC 12) . 
1FMUC12) *WC12) ) 

NRAY-NRAY-l 
60 TO 202 

c 

C EXIT 

C 

990 CALL DVCHK1JJ) 

GO TO 1992, 995 >.JJ 

992 WRITE (6.9931 

993 FORMAT 1 34 HOD I V I OE CHECK ON WHEN EXIT REACHCD/1H1) 
CALL EXIT 

995 WRITE 16.996) 

996 FORMAT! I 2H0EN0 OF CASE/1 Hi) 

997 LINE-0 
CALL PREND 
CALL EXIT 

C 

C ERROR 

9 

998 WRITE <6*999) MODE 

999 FORMAT 122HOPROGRAM ERROR IN M0DE.t2/im* 

GO TO 997 

END • 


SUBROUT I NE GENL ( GXA . GY A . GT A , GMUA , GW A . GXB . GyB ♦ GTB • 

1 GMUB • GWB . GXC . GYC . GTC . GMUC . GWC ) 

C 

C GENERAL POINT SUBROUTINE 

C 

DIMENSION CASE 1 14 ),XA ( 400 ).YA (400 l.THETA 1400). FMUA (400). 
1 WA (400) .CHARAC (6.700 ) ,XC (400 ) .YC1400 ) .THETC 1400) * 

2FMUC (400 ) .WC1400 ) .FMACHC (400 ) 

COMMON CASE. XA.YA. THETA. FMUA.WA. 

1 FMACHC . XC . YC . THETC .FMUC . WC . CHARAC * 

2GM 1 , GM 1 2 . XP . L I NE . NA . NC , NRA Y . NREFL 
AS1N!X)-ATAN2(X.S0RT(1 ,-X**2 )) 

SI -5 IN (GMUB) 

B2-COS 1 GTB— GMUB ) 

B3-5 I N{ GMUA ) 

B4-C0S (GTA+GMUA ) 

1 D1-SIN(GTA+GMUA»/B4 
D2-S IN ( GTB-GMUB )/B2 
D4-03/COS ( GMUA ) 

05-81 /COS (GMUB) 

C EQUATION (C3 ) 

GXC- (GXA*DI -GYA+GYB-GXBWD2 )/ <DI -02 ) 

D3-GXC— GXA 
D6-D3/GYA 

C EQUATION 1C4) 

8 GYC-D3*D1+GYA 
C EQUATION (CIO) 

FMB- (81 *SINCGTB )*D5 ) /B2 

C EQUATION ( C9 ) 

FLA-(B3*SIN(GTA)*D4)/B4 
C I — GTA*D4+D6*FLA 
C2-GTB#05+ 1 GXC-GXB ) 4FM8/GYB 
C EQUATION (CB) 

GTC* ( -GWA-GWA4C 1 +GWB+GWB*C2 ) / ( 6WA*D4+GWB*05 ) 

C EQUATION (CU) 

GWC»GWA+GWA* (D4*(GTC-GTA )+06*FLA > 

Cl -GM1 2* ( 1 . 0/GWC**2- 1*0) 

C EQUATION (C12) 

GMUC- AS INI SORT (Cl ) ) 

28 RETURN 
END 



SUBROUT I NE CENTL ( GXA .6YA . GTA .GMUA . GWA ♦ 
lGXC.GMUC.GWC) 

C 

c centerline: point subroutine 

c 

DIMENSION CASE (14 ) tXA (400 ) « YA (*00 ) . THETA ( 400 > « FMUA (400 > • 
1 WA (400 ) . CHARAC ( 6 • 700 ) « XC ( 4 00 ) ♦ YC ( A 00 ) « THETC < AOO » . 

2FMUC (400 ) iHC ( AOO > .FMACHC 14 00 > 

COMMON CASE.XA. YA, THETA. FMUA ,WA. 

1 FMACHC . XC . YC ♦ THETC , FMUC . WC , CHARAC . 

2GM1 .GM12.XP.LINE.NA » NC «NRA Y . NREFL 
ASIN(X)«ATAN2(X.SQRT(1 «-X**2) ) 

BI »SIN(GMUA ) 

B2*C0S ( GTA+GMUA > 

1 DI «S!N(GTA+GMUA >/B2 
02-B1/C0S ( GMUA ) 

C EQUATION (C28) 

GXC-GXA— GYA/D1 
C EQUATION (C9> 

FLA- (B I *S INI GT A >#D2 )/B2 
C EQUATION (CM ) 

GWC-GWA+GWA* (-GTA*D2+ C GXC-GXA ) +FLA/GYA ) 

CI-GM12*U • 0/GWC**2— 1 .0) 

C EQUATION (Cl 2) 

GMUC-AS IN (SORT (Cl ) ) 

18 RETURN 
END 


SUBROUTINE OFCNT ( GXA »GYA »GTA .GMUA « GWA »GXB , 

! GMUB . GWB . GXC , GYC.GTC »GMUC * GWC ) 

C 

C 1 POINT OFF CENTERLINE SUBROUTINE 

C 

DIMENSION CASE (14 J , XAMOO) .YA(400 J .THETA (400 ) .FMUA (400) . 
1 WA (400) .CHARAC (6 . 700 ) .XC (400 ) *YC (400 }. THETC ( 400 ) « 

2FMUC (400) » WC ( 400 ) .FMACHC (400) 

COMMON CASE. XA.YA. THETA. FMUA .WA, 

1 FMACHC t XC . YC .THETC .FMUC , WC . CHARAC . 

2GMI ,GMI 2.XP.LINE.NA.NC.NRAY.NREFL 
AS I N (X ) - AT AN 2 (X. SORT ( 1 .-X**2 > ) 

Bl -S IN ( GMUA ) 

B2-COS (GTA+GMUA ) 

I 01 =B1 /COS (GMUA ) 

D2-S I N ( GMUB ) /COS ( GMUB ) 

D4-SIN (GTA+GMUA )/B2 
DS-S I N { —GMUB ) /COS (-GMUB ) 

C EQUATION (C3 ) 

GXC- (GXA + D4-GYA— GXB*D5 >/(D4-D5 ) 

D3-GKC-GXA 

C EQUATION (C4) 

GYC-D3*D4+GYA 

C EQUATION (C9> 

FL*» (B 1 *S I N ( GTA )*D1 )/B2 
Cl *— GTA#D I +D3+FLA/GYA 
C EQUATION (C27) « 

GTC- ( — GWA — GWA*C 1 +GWB )/ ( GWA+D 1 +2 * *GWB*D2 ) 

C EQUATION (Cll) 

GWC-GWA+GWA* (D 1 * ( GTC-GTA ) +D3*FLA/GYA ) 

Cl -GMI2*( I • 0/GWC**2— 1 .0) 

C EQUATION (C 1 2 ) 

1 4 GMUC-AS I N ( SORT { Cl ) ) 

18 RETURN 
END 


SUBROUTINE BNORY (BXA.BYA.BTA.BMUA.BWA.BXB.BYB.BTB.BMUB. 
1BWB.BXC.BYC. BTC. BMUC .BWC) 

C 

C BOUNDARY POINT SUBROUTINE 

C 

B 1 -COS ( BTB-BMUB ) 

B2-S IN ( BMUB ) 

1 Dl *SIN(BTA J/COS (BTA ) 

D2*S I N < STB -B MUB ) /B 1 
D3*S I N ( BMUB ) /COS ( BMUB ) 

C EQUATION (CIB) 

BXC- (BXA*DI-BYA+BYB-BXB*D2 )/ (D1-D2 ) 

C EQUATION (Cl 9 ) 

BYC- ( BXC— BXA ) *D1 +BYA 

BMUC-BMUA 

BWC-BWA 

C EQUATION (CIO) 

FMB- ( B2*S I N ( BTB )*03)/Bl 
C EQUATION (C20 ) 

BTC-BTB+ ( — (BWC-BWB )/BWB+( (BXC-BXB ) /BYB ) *FM8 ) /D3 
12 RETURN 
END 


SUBROUTINE SAMFM f SXB . SYB , STB • SMUB. SWB . 

1 SXA, 5YA . STA . SMUA , SWA .SXC.SYC . STC , SMUC . SWC > 

C 

C SAME FAMILY POINT SUBROUTINE 

C 

DIMENSION CASE ( 14 > .XA(400 » .YA (400 ) .THETA (400 ) .FMUA (400 > . 
1 WA (400 ) .CHARAC (6.700 )«XC(400 ).YC(400 >. THETC ( 400 ) . 

2FMUC ( 400 ) « WC ( 400 ) .FMACHC (400 ) 

COMMON CASE. XA.YA. THETA. FMUA «WA. 



1 FMACHC . XC « YC . THETC . FMUC . WC . CHARAC » 

2GMI «GM) 2 • XP.L I NE »NA < NC *NPA Y * NREFl 
AS I N (X ) -ATAN2 (X.SQRT f TTT > > 

Bl -SIN(SMUB) 

B2-COSCSTB+SMUB) 

B3*S!N(SMUA> 

B4«C0S(STA+SMUA ) 

01 -S I N ( ST A+SMUA >/B4 
D2-S ! N ( STB+SMUB >/B2 
D4«B3/C0S(SMUA ) 

D5*81/C0S(SMUB> 

C EOUATION (C23 ) 

SXC*(5XA*DI-SYA+SYB-SXB*D2 )/(Dl-D2) 

D3-SXC-SXA 

C EQUATION <C4 ) 

SYC-D34DI+SYA 

C EQUATION (C24I 

FMB*(BI*SIN(STB)*D5 1/82 
C EQUATION (C9> 

FLA»(B3*SIN(STA>*04)/B4 
Cl *-STA*D4+D3*FLA/SVA 
C2-STB*D5+f CSXC-SXB)*FMB)/SYB 
C EQUATION fCB) 

STC-(-SWA-SWA*Cl+SWB + SWB*C2 >/<SWA*D4+SWB*D5> 

C EQUATION (Cl!) 

SWC-SWA+SjfA* (D4*(STC-STA I+D3+FLA/SYA > 

Cl *GM!2«(! .0/SwC**2-l .0) 

IF (C 1 >97.24.24 
24 TT-SQRTCC1) 

TTT" I •— TT**2 
I F ( TTT > 97 . 26 * 26 
C EQUATION (CI2) 

26 SMUC-ASIN(TT) 

GO TO 28 

97 WRITE <6. 98) 

98 FORMAT (37H SORT ERROR IN SAMFM POINT SUBROUTINE) 
28 RETURN 

END 


SUBROUTINE TEST 

C 

C SUBROUTINE TO TEST FOP CROSSINGS 

C USES 3 SAME family POINTS 

C 

DIMENSION CASE (14).XA( 400 )*YA(400). THETA (400 ) ,FMUA(400) . 

1 WA (400 ) .CHAR AC (6. 700 ) . XC ( 400 ) . YC (400 ) • THETC ( 400 > . 

2FMUC (400 ) ,WC (400 > .FMACHC (400 > 

COMMON CASE. XA.YA. THETA. FMUA.WA. 

1 FMACHC • XC . YC .THETC , FMgC . WC i CHARAC « 

26M 1 *GMI2,i<P.L!NE.NA*NC.NRAY,NReEL 
500 CONV-57.29578 

CON* THETC (NC— I )-FMUC (NC-i ) 

SLOPE »S I N { CON ) /COS C CON ) 

CB*YC(NC-1 )-XC (NC-I >*SLOPE 
XASAV-XA (NA+1 ) 

NASAV*NA 
GO TO 515 
C 

C COMPUTE NEW GENL 

C 

510 CALL GENL (XA (NA ), YA ( NA )* THETA (NA > .FMUA (NA).WA(NA). 

1 XC (NC— 1 >,YCfNC-l ) .THETC (NC-1 > ,FMUC(NC-1 ) *WC(NC-1 ) • 

2XC (NC). YC(NC) . THETC (NC ) «FMUC (NC ) « WC ( NC ) ) 

c 

C COMPUTE 1 PT 

C 

515 CALC SAMFM (XA (NA > .YA (NA) .THETA (NA ) .FMUA (NA) ,WA (NA) . 

1XA (NA+1 ) .YA (NA+1 ) *THET A (NA + 1 ) . FMUA (NA+1 ) . WA ( NA+1 ) . 

2X1 • Y1 .THET1 « FMU 1 .Wt ) 

c 

C COMPUTE 2 PT 

C 

520 CALL SAMFM(XA (NA+1 > ,YA (NA+J ) .THETA (NA+ 1 ) .FMUA (NA+1 ) *WA(NA+1 ) 
1 XA (NA+2 ) ,YA (NA+2 ) .THETA (NA+2 ) «FMUA(NA+2) *WA(NA+2 ) . 

2X2 . Y2 • THET2 . FMU2 . W2 ) 

C 

C COMPUTE 3 PT 

c 

CALL SAMFM (XA (NA+2 ) , YA (NA+2 > .THETA (NA+2) .FMUA (NA+2) .WA (NA+2) 
1 XA (NA+3 ) .YA(NA+3) . THETA (NA+3 ). FMUA (NA+3) .WA(NA+3) . 

2X3. Y3 . THET3 . FMU3 . W3 ) 

C 

C FtNO X ANO Y MIN ANO MAX FC** BOX 

C 

525 XMIN*AMIN1 (XCCNC) .XA (NA) «XA(NA+1 ) * XA (NA+2) ,XA (NA+3)) 

526 XMAX-AMAX1 (XC(NC) ,XA (NA) tXA(NA+l > .XA (NA+2 ). XA ( NA+3 ) > 

527 YMIN-AMIN1 (YC(NC) , YA (NA) . YA (NA+I ) .YA (NA+2) . YA (NA+3) ) 

528 YMAX-AMAX1 (YC(NC) ,YA (NA) ,YA(NA+1 ) .YA (NA+2 ) .YA (NA+3 ) ) 

C 

C CHECK t PT FOR CROSSING 

C 

call crosscxi *yi .xmax.xmin.ymax*ymin. 

I THETA (NA) .FMUA ( NA ), THETA (NA+ 1 ) .FMUA ( NA+1 ) » 

2SLOPE.C8. IT) 

IT- IT 

GO T0f530.535) . IT 
530 KOOE-2 

GO TO 538 
535 KOOE-1 

6 



C CHECK 2 PT EOR CROSSING 

C 

538 CALL CROSS (X2« Y2*XMAX»XM IN * YMAX • YM INi 

I TMET A ( NA+ 1 >. FMUA (NA+1 ) .THETA (NA+2 > «FMUA (NA+2 > » 
2SLOPE.CB. IT) 

IT-1T 

GO TO ( 540 • 545 ) « IT 
540 GO TO(542.544) .KODE 
542 KODE-3 

GO TO 545 

544 KODE-4 
C 

C CHECK 3 PT FOR CROSSING 

C 

545 CALL CR0SSCX3.Y3.XMAX*XMIN.YMAX« YMIN. 

1 THETA t NA+2 ) « FMUA (NA+2 ) « THETA (NA+3 ) * FMUA (NA+3 ) * 
2SLOPE«CB. IT) 

IT-tT 

GO TO ( 550 *560 ) « IT 
C 

C FIND COMBINATION OF CROSSINGS TO CHECK 

C 

550 GO TOf 750*760*770*780) *KOOE 
500 GO TO 171 0*720 *730 *740 ) »K0DE 
C 

C- NO CROSSING 

C 

710 L INE-L INE— (NA— NASAV ) 

RETURN 

C 

C I PT IS CROSSING TO USE 

C 

720 TPR-THETl*CONV 
FMUPR-FMU1 4C0NV 

WRITE (6.7221X1 «Y1 .TPR.FMUPR.W1 
722 FORMAT (6H 4* I ** I 1 X* 3EJ 7. 8 « 2E 1 6 *8 ) 

IF (NC-NRAY 1725.725* 724 

724 NREFL-NREFL-I 
GO TO 726 

725 NRAY-NRAY-I 

726 NA-NA+1 
XA (NA ) - XI 
YA(NA)*Y1 
THETA (NA)-THETl 
FMUA(NA)-FMU1 
WA(NA)*W1 

GO TO 310 
C 

C 2 PT IS CROSSING TO USE 

C 

T*9 TPQ«TH*T**<$©NV 

FMUPR«^MU2*C0NV 

WRITE <6*732 1X2. Y2*TPR.FMUPR*wa 
732 FORMAT (6H *42*4 11 X. 3£I 7. 8 * 2E 1 6.8 ) 

IF (NC-NRAY >735*735.734 

734 NREFL-NREFL-1 
GO TO 736 

735 NRAY-MRAY-1 

736 XAtNA+1 >-XA(NA ) 

YA (NA+1 ) -YA ( NA ) 

THETA (NA+1 ) -THETA (NA > 

FMUA (NA+1 l-FMUA (NA) 

WA (NA+1 )>WA(NA ) 

NA-NA+1 

XA (NA+1 )»X2 
YA (NA+1 ) -Y2 
THETA < NA+1 ) — THET2 
FMUA (NA+1 > -FMU2 
WA (NA+1 )-W2 
GO TO 515 
C 

C BOTH 1 AND 2 ARE CROSSINGS* SELECT ONE TO USE 

C 

740 DIFFl-ABSfXl-XASAV) 

0IFF2-AB5 (X2— XASAV ) 

IF (0IFF1 — DIFF2 >720.720*730 
C 

C 3 PT IS CROSSING TO USE 

C 

750 TPR-THET34C0NV 
FMUPR-FMU34C0NV 

WRITE (6.732)X3*Y3.TPR,FMUPR.W3 
752 FORMAT ( 6H 4*3** 1 1 X « 3E1 7. 8, 2E1 6 .8 > 

IF (NC-NRA Y > 755.755 . 754 

754 NREFL-NREFL-1 
GO TO 756 

755 NRAY-NRAY-1 

756 XA (NA+2 >-XA (NA+1 ) 

YA (NA+2 >*YA (NA+1 ) 

THETA (NA+2) -THETA (NA+ I ) 

FMUA (NA+2 ) -FMUA (NA+1 ) 

WA (NA+2 )— WA (NA+1 ) 

XA (NA+1 )-XA ( NA > 

YA (NA+1 )-YA (NA ) 

THETA (NA+1 ) -THETA (NA ) 

FMUA ( NA+ 1 ) -FMUA ( NA > 

WA (NA+1 )«WA(NA ) 

NA-NA+1 
XA (NA+2 ) -X3 
YA (NA+2 ) — Y3 
THETA (NA+2 1-THET3 
FMWA(N4+2|-FMW3 





WA (NA + 2 )«W3 
c GO TO 820 

C BOTH 1 AND 3 ARE CROSSINGS* SELECT ONE TO USE 

C 

760 DIFE1 -ABS (XJ —XASAV) 

DIFF3-ABS (X3-XASAV) 

IF CDIFF1 —0 1 FF3 ) 720 « 720 • 750 
C 

C BOTH 2 AND 3 ARE CROSSINGS. SELECT ONE TO USE 

C 

770 D I FF2- ABS ( X2— XASAV ) 

0IFF3-ABS ( X3— XASAV ) 

IFtDIFF2-DIFF3l730.730.750 

C 

C 1.2. and 3 ARE CROSSINGS. SELECT ONE TO USE 

C 

780 DIFFJ -ABS(Xt-XASAV) 

O I FF2= ABS ( X2— XASAV ) 

D IFF3* ABS ( X3— XASAV ) 

IF (DIFF1-DIFF2 >782,782.785 

782 I F f D ! FF 1 — D I FF3 ) 720 . 720 . 750 

785 IF (DIFF2-DIFF3) 730. 730.750 
END 


SUBROUTINE TEST2MSIG) 

C 

C MERGES 2 CHARACTERISTIC LINES IF THEY CROSS 

C I S I G* 1 FOR NO CROSSING 

C ISIG-2 FOR CROSSING FOUND 

C 

D I MENS ION CASE ( I 4 I « XA ( 400 ) . Y A (400 ) . THETA 1 400 ) . FHUA ( 400 ) . 
1 WA (400 ) .CHARAC (6.700 ) . XC (4 00 ) . YC ( 4 00 ) . THETC I 400 > . 

2=-MUC (400) .WC(400) .FMACHC (400 ) 

COMMON CASE. XA.YA, THETA. FMUA »WA« 

I FMACHC . XC . YC . THETC . FMUC . WC . CHARAC . 

2GM1 ,GMI 2.XP.LINE.NA .NC .NR AY . NREFL 
ANGLE = — 1 .5^07063 

1 IF (THETA (NA I+FMuA (NA J-ANGLF ) 3.2.2 

2 IF (XC(NC)-XA (NA ) > 6.5.5 

3 IF ( XC ( NC ) — XA (NA ) ) 5.5.6 
C 

C NO CROSSING 

C 

5 IS I G* 1 
RETURN 

C 

C CROSSING FOUND 

C 

6 CALL SL1TET(2.KOOOFX) 

GO TO (8. I 0) .KOOOFX 
8 CALL SLITE (2) 

GO TO 15 
10 NREFL -NREFL- 1 
C 

C MOVE REMAINING A ARRAY TO C ARRAY AND PRINT 

C 

15 NA-NA+1 
NC-NC+l 

20 XCINC-I >*XA (NA-1 > 

YC (NC— 1 }-VA(NA-l ) 

THETC (NC-1 ) -THETA (NA-1 ) 

FMUC(NC-1 ) — FMUA (NA— 1 ) 

WC (NC-1 ) -WA (NA-1 ) 

FMACHC (NC-1 >■! ./S IN (FMUC (NC-1 ) ) 

CALL PRINT (FMACHC (NC-1 I.XC(NC-1 >.YC(NC-1 ).THETC(NC-1 ) * 
l FMUC (NC-1 ) .WC (NC-1 ) ) 

L1FOLINE- I 
IF (LINE! 998.25.15 
25 ISIG-2 
NC-NC-1 
RETURN 

998 WRITE (6.999) 

999 FORMAT (20H0ERR0R IN TEST2 SU8R/1H1) 

CALL EXIT 

END 


SUBROUTINE CROSS (X.Y.XMAX. XM IN. YMAX . YM IN. 

I T1 , FM I . T2 , FM2 , SLOPE . CB , I T ) 

C 

C CHECKS A PARTICULAR POINT TO SEE IF A CROSSING 

C IT-1. CROSSING FOUND 

C IT-2. CROSSING NOT FOUND 

C 

C CHECK b-IMIT BOX 

C 

IF (Y-YMINI25.25.2 
2 I F ( Y— YMAX ) 4 « 25 ♦ 25 
4 I F ( X— XM IN) 25 * 25 *6 
6 IFtX— XMAXJ8.25.25 
C 

C CHECK CONVERGING LINES 

C 

8 IF t (T1 +FM1 )- (T2+FM2 ) ) 10. 10.25 
C 

C CHECK IF PAST LINE 8 


10 CSF-Y— SLOPE* X 



12 lF(SL0PEM4ti6*IA 
t 

C NEGATIVE SLOPE 

C 

14 IF (CSF—CB ) 18. 1 8 <25 

C 

C POSITIVE SLOPE 

c 

16 IF (CSF—CB ) 25 *18*16 
C 

C CROSSING FOUND 

C 

18 I T— 1 
RETURN 
C 

C NO CROSSING FOUND 

C 

29 IT-2 
RETURN 
END 


SUBROUTINE PRINT (PM .PX .PY. PT .PMU.PW ) 

C 

C PRINTS A POINT IN THE NETWORK 

C 

DIMENSION CASE C 1 4 ) • XA (400 ) • YA (400 ) .THET A (400 > »FMUA (400 ) . 
I WA (400). CHARAC (6. 700 ).XC(400) *YC( 400) .THETC <4001 . 

2FMUC C400 > .WC (400 ) .FMACHC ( 400 I 
COMMON CASE. XA.YA. THETA. FMUA.WA. 

I FMACHC • XC . YC . THETC . FMUC • WC . CHARAC . 

2GM1 .GM12.XP.LINE.NA.NC.NRAY.NREFL 
DEG-97. 29978 
TEMP-1 .+GM12*PM**2 
C EQUATION (C 1 4 ) 

TTZ-1 ./TEMP 

C EQUATION (C 1 9 ) 

RRZ-1 .ZTEMP»*XP 
TPR-PT *DEG 
FMUPR- PMU*DEG 

WRITE(6.10)PM,PX.PY.TPR.FMUPR.PW.TTZ.RRZ 
10 FORMAT ( 4EI 7*8 .4EI 6*8) 

RETURN 

END 


SUBROUTINE PREND 
C 

C MARKS ENO OF A CHARACTERISTIC LINE AND 

C WfclTES A TA**E FOR OPTIONAL PLOYtthlA 

C 

DIMENSION CASE (14). XA (400 ) . YA (400 ) . THETA ( 400 ) . FMUA ( 400 > » 
1 WA (400 ). CHARAC (6, 700 ). XC(400 I. YC (400) .THETC (400). 

2FMUC (400 ) .WC 1 400 ) .FMACHC ( 400 ) 

COMMON CASE. XA.YA. THETA. FMUA.WA. 

I FMACHC . XC . YC . THETC . FMUC . WC . CHARAC . 

2GM1 .GMJ2.XP.LINE.NA.NC.NRAY.NREFL 
IF (LINE )2 . 2 • 4 
2 WRITE ( 9 >L INE 
END FILE 9 
RETURN 

4 WRITE(6.B)NRAY.NREFL 
8 FORMAT ( 1H0.2 I6/'l HI ) 

C 

C WRITE PLOTTING TAPE 

C 

WRITE ( 9 )L INE 

WRITE(9> (FMACHC 1 1 ).XC( I ) *YC( I ) • 1 -1 .LINE) 

RETURN 

END 


ENTRY .UN09. 

-UN09. PZE UN1T09 

UNIT09 FILE »B ( 1 > .OUTPUT. BLOCK-296. B IN. HOLD 

END 
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Schlieren mirror 



(a) Test setup in 4l-foot vacuum sphere. 
Figure 1.- Schematics of test facilities. 


52 






..250 






jet boundary 




56 


Figure 3*- General description of exhaust jet. 



FACTORS AFFECTING ROCKET EXHAUST SIZE AND SHAPE 


M. Nozzle-exit Mach number 

J 

Nozzle half- angle 

7 Ratio of epedfic heats 

Pj Static pressure at nozzle exit 

p^ Ambient static pressure 

“l, - V 1 - V N - 9 N 


where 

Prandtl-Meyer expansion angle for boundary Mach number 
Prandtl-Meyer expansion single for nozzle-exit Mach number 



Figure 4.- Method of calculating initial turning angle of jet boundary at the nozzle lip. 








Figure 7.- Comparison of jet boundary coordinates for foldback procedure with those for modified method in which foldback is 

eliminated. Mi = 1 . 0 . 






(a) Expanded scale showing boundaries near nozzle exit for 




(b) Complete boundaries 







(c) Expanded scale showing boundaries near nozzle exit for 










926 Concluded 



Figure 8.- Continued. 





(g) Jet boundaries and internal shocks for 




Figure 9 *- Location of maximum jet exhaust boundary coordinates in terms of nozzle radius 







Figure 9.- Concluded. 
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Figure 10.- Concluded 
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Figure 11 .- Continued. L-6^-395 








83 


L-64-398 
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Figure 11.- Continued. L-64-400 





Figure 11.- Continued. L-64-402 
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Figure 11.- Continued. L-64-4o4 






Figure 12.- Schlieren photographs of jet exhaust plumes with superimposed, theoretical boundaries and internal shock. 
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Figure 12.- Concluded. L-64-407 
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Figure l4.- Comparison of "boundaries computed "by characteristic method with those computed by empirical method. 





Figure 15.- Effect of jet pressure ratio upon nondimensional distance along jet axis from plane of nozzle exit to location 

of Riemann wave. 







Figure 15.- Continued 






Figure 15.- Concluded, 




